# Multilateration and Signal Matching with Unknown Emission Times

Mireille Boutin  Department of Mathematics, Purdue University, 150 N.
University St., West Lafayette, IN, USA 47907<EMAIL_ADDRESS>and  Gregor
Kemper  Technische Universität München, Zentrum Mathematik - M11,
Boltzmannstr. 3, 85748 Garching, Germany<EMAIL_ADDRESS>
(Date: June, 2022)

###### Abstract.

Assume that a source emits a signal in $3$-dimensional space at an unknown
time, which is received by at least $5$ sensors. In almost all cases the
emission time and source position can be worked out uniquely from the
knowledge of the times when the sensors receive the signal. The task to do so
is the multilateration problem. But when there are several emission events
originating from several sources, the received signals must first be matched
in order to find the emission times and source positions. In this paper, we
propose to use algebraic relations between reception times to achieve this
matching. A special case occurs when the signals are actually echoes from a
single emission event. In this case, solving the signal matching problem
allows one to reconstruct the positions of the reflecting walls. We show that,
no matter where the walls are situated, our matching algorithm works correctly
for almost all positions of the sensors.

In the first section of this paper we consider the multilateration problem,
which is equivalent to the GPS-problem, and give a simple algebraic solution
that applies in all dimensions.

###### Key words and phrases:

GPS-problem, multilateration, echo sorting, TDOA disambiguation, shape
reconstruction

###### 2010 Mathematics Subject Classification:

51K99, 13P10, 13P25

## Introduction

Positioning is a ubiquitous problem in engineering. For example, one may want
to determine the location of an object such as a vehicle, locate an event such
as an earthquake, calibrate an array of devices such as microphones, or draw
the map of an environment such as a building. In many scenarios, the objects
to be located can emit a signal. In such case, one can use an array of
receivers with known geometry to determine the objects location with respect
to the position of the receiving array. Alternatively, the objects may be
equipped with a receiver so to be located with an array of sources with known
geometry. While both problems are dual to each other, their difficulty and
conditioning can vary significantly depending on the specific setup scenario
and constraints imposed. In some cases, the underlying mathematical problem
may actually be ill-posed.

This paper is concerned with multilateration, which is the task of determining
the position of one or more sources emitting a wave signal (e.g.,
electromagnetic, acoustic, or seismic waves). More specifically, we are trying
to determine the position of sources sending out a signal from measurements of
the times when this signal is received by various sensors situated at known
positions. We assume that the clocks on the receivers are synchronized
together, but not with the clock of the sources. In other words, the time of
signal emission is unknown to the sensors, and thus the differences of
arrivals (TDOAs) are the only meaningful available information. Therefore one
also speaks of pseudo-range multilateration.

The literature on the well-posedness of multilateration problems is sparse. As
far as we know, even the well-known GPS positioning problem, which we analyze
in Section 1, has not been thoroughly studied. This corresponds to the problem
of determining the position of one emission event (at an unknown time)
received by several (synchronized) sensors at known locations. It turns out
that, even in this simple scenario, the TDOAs may not uniquely determine the
position of the source.

More generally, we consider the case of several sources emitting
undistinguishable signals at unknown times (e.g., earthquakes
[?sha2021reliability] or gunshots). Clearly the order in which the signals
arrive at the sensors can be wildly different from the order in which the
signals were emitted. Therefore before feeding the reception times into a
multilateration algorithm, a matching must be performed, identifying those
received signals that come from the same emission event with each other.
Ideally, this matching process should also discard spurious signals that are
registered by just one sensor, i.e., not match such signals to any others.
Note that a mistake in the matching process, or accidentally including a
spurious signal in a match, will result in determining an emission event (time
and position of emission) which never took place.

The problem of matching the sound events produced at a known time was
previously studied in [?DPWLV1, ?jager2016room,
?el2017time,?Boutin:Kemper:2019]. In the case where the times of emission are
unknown, the problem is known as “TDOA disambiguation.” Two sources of TDOA
ambiguity are considered in the literature. The first one is the “multipath
ambiguity” which is caused by the reverberation of the signal on objects in
the environment and leads to spurious events. The triangle inequality, a zero
cyclic sum condition, and characteristics of the cross-correlation and
autocorrelation of the microphone signals are used to help disambiguate such
cases in [?kreissig2013fast,?zannini2010improved]. The second one is the
“multisource ambiguity” resulting from multiple sources emitting the same
signal. A few of the false matches due to this can be ruled out using the
triangle inequality. But, as far as we know, a more rigorous criterion for the
case where the times of emission are unknown has not previously been proposed.

In this paper we show that if there are (at least) $5$ sensors in
$3$-dimensional space, then the reception times (more precisely, TDOAs) of a
signal coming from a single emission event satisfy a certain algebraic
relation. We propose to use this relation to perform the signal matching. So
if a selection of signal reception times, one for each sensor, satisfies this
relation, then these reception times are accepted as coming from the same
event (see Algorithm 2.3). Moreover, spurious signals registered by just one
sensor will almost certainly not satisfy the relation and will therefore not
be included in a match.

Of course the ”almost” in the previous sentence is an important issue. In
fact, there is no way to rule out the possibility that a spurious signal is
registered at such an unlucky time that our matching algorithm, or any other
algorithm based on the available information, falsely includes this signal in
a match. This is true not only for spurious signals but also for signals
coming from a real emission event that happened at an unlucky time. In
particular, there is no way to position the sensors such that this possibility
can be ruled out.

The situation becomes different, however, if the emission events are in fact
just echoes from a single event. More precisely, consider an arrangement of
flat surfaces (”walls”) that reflect a sound signal emitted from an
omnidirectional loudspeaker. Assuming the signal is of high frequency, we use
the ray acoustics approximation. This means that the signals are virtually
emitted from the point given by reflecting the loudspeaker position at the
walls, and all (virtual) emissions happen simultaneously. Now we bring in $5$
microphones at known positions. These record the echoes of the sound emission
and feed them into the matching algorithm. After that, the wall positions can
be determined by multilateration. Notice that in contrast to our previous
paper [?Boutin:Kemper:2019], we do not assume that the loudspeaker and the
microphones have synchronized clocks and communicate the times of signal
emission time. So the common emission time is still unknown, and pseudo-range
multilateration is required. In Theorem 3.1 in this paper, we show that in
this situation almost all microphone positions are good, in the sense that no
false matches can happen. As explained above, this is in contrast to the
situation where the (virtual) emission events are not assumed to be
simultaneous. The proof of the theorem uses methods from computational
commutative algebra. This is something that it has in common with the proof of
the main results from [?Boutin:Kemper:2019]. However, when we designed the
proof of Theorem 3.1 we were surprised to find that the difficulties that
arose were quite different from those in [?Boutin:Kemper:2019].

Even though the emphasis of this paper may lie on the matching problem and on
multiple wall detection, we also study the pseudo-range multilateration
problem itself. In fact, we present a simple algebraic solution algorithm. For
simplicity, we formulate this in three dimensions, but it really works for all
dimensions $>1$. We give a self-contained and rigid proof for the validity. It
is well-known that with just $4$ sensors, the pseudo-range multilateration
problem usually has $2$ solutions. We show by examples that this may also
happen if there are $5$ sensors, even if no $4$ of them are coplanar. This is
not a shortcoming of our algorithm: in our example, the available information
of TDOAs simply does not allow to disambiguate the solutions.

The paper is organized in three sections. The first section introduces the
notation and discusses the pseudo-range multilateration problem, giving the
solution algorithm in Theorem 1.1. In the second section we then turn our
attention to the case of multiple emission events. We present and prove the
relation that holds between reception times coming from the same event, and
derive the matching algorithm (Algorithm 2.3) from this. The final section
deals with the situation of matching echoes from a single sound event. The
main result (Theorem 3.1) from that section says that, loosely speaking,
almost all microphone positions are good.

Acknowledgments. A great part of this work was done during a research stay of
the authors at the Mathematical Sciences Research Institute (MSRI) in Berkeley
within the 2022 Summer Research in Mathematics program. We thank the MSRI team
for creating a uniquely stimulating atmosphere at the institute. The program
provided us a with a perfect research environment and a chance to concentrate
fully on getting things done. This work has also benefited from a research
stay of the two authors at the Banff International Research Station for
Mathematical Innovation and Discovery (BIRS) under the “Research in Teams”
program. We would like to thank BIRS for its hospitality and for providing an
optimal working environment. We also thank Stefan Wetkge and Timm Oertel for
stimulating discussions.

##  1\. Multilateration and the GPS-problem

In this section we look at the pseudo-range multilateration problem, which is
equivalent to the GPS-problem (see below), and present two simple algebraic
(“direct”, as opposed to iterative) solutions, see Theorem 1.1. This topic has
received considerable interest in the literature (see ?Bancroft:1985
[?Bancroft:1985], ?Krause:1987 [?Krause:1987], ?Chaffee:Abel:1994
[?Chaffee:Abel:1994], ?li2010design [?li2010design], ?lundberg2001alternative
[?lundberg2001alternative], and ?Beck:Pan:2012 [?Beck:Pan:2012]). But our
results are general and appear to be new.

Our primary interest is in the following situation: A source at an unknown
position $\mathbf{x}\in{\mathbb{R}}^{3}$ emits a signal, in our applications
usually by sound, at an unknown time $t$. (In fact, everything we are about to
say can easily be adapted to ${\mathbb{R}}^{n}$ with $n\geq 2$, but not to
${\mathbb{R}}^{1}$, see Remark 1.2.) There are $m$ sensors at known positions
$\mathbf{a}_{1},\ldots,\mathbf{a}_{m}\in{\mathbb{R}}^{3}$. They receive the
signal at times $t_{1},\ldots,t_{m}$. We choose the unit of time such that the
signal propagation speed becomes $1$. So we have

$\lVert\mathbf{a}_{i}-\mathbf{x}\rVert=t_{i}-t\qquad(i=1,\ldots,m)$ (1.1)

The task now is to work out the position $\mathbf{x}$ and the emission time
$t$. The very same equations arise if there are $m$ sources at known positions
$\mathbf{a}_{i}$ emitting signals at known times, which are then received by a
device at an unknown position $\mathbf{x}$. In this case the $t_{i}$ are the
differences between the reception times according to the clock on the device
and the emission times according to the (near-perfect) clocks on the sources,
and $t$ is the (unknown) bias between the clock on the receiver and the clocks
on the sources. This is the GPS-problem.

We will work with the slightly weaker equations

$\lVert\mathbf{a}_{i}-\mathbf{x}\rVert=|t_{i}-t|\qquad(i=1,\ldots,m).$ (1.2)

Writing $L:=\left(\begin{smallmatrix}-1\\\ &1\\\ &&1\\\
&&&1\end{smallmatrix}\right)$,
$\widetilde{\mathbf{a}_{i}}:=\left(\begin{smallmatrix}t_{i}\\\
\mathbf{a}_{i}\end{smallmatrix}\right)$ and
$\widetilde{\mathbf{x}}:=\left(\begin{smallmatrix}t\\\
\mathbf{x}\end{smallmatrix}\right)$, we have

$\lVert\mathbf{a}_{i}-\mathbf{x}\rVert^{2}-(t_{i}-t)^{2}=(\widetilde{\mathbf{a}_{i}}-\widetilde{\mathbf{x}})^{T}\cdot
L\cdot(\widetilde{\mathbf{a}_{i}}-\widetilde{\mathbf{x}})=\widetilde{\mathbf{a}_{i}}^{T}L\widetilde{\mathbf{a}_{i}}-2\widetilde{\mathbf{a}_{i}}^{T}L\widetilde{\mathbf{x}}+\widetilde{\mathbf{x}}^{T}L\widetilde{\mathbf{x}}=\\\
\lVert\mathbf{a}_{i}\rVert^{2}-t_{i}^{2}+2t_{i}t-2\mathbf{a}_{i}^{T}\mathbf{x}+\lVert\mathbf{x}\rVert^{2}-t^{2},$

so (1.2) is equivalent to

$-2t_{i}t+2\mathbf{a}_{i}^{T}\mathbf{x}-\lVert\mathbf{x}\rVert^{2}+t^{2}=\lVert\mathbf{a}_{i}\rVert^{2}-t_{i}^{2}\qquad(i=1,\ldots,m).$
(1.3)

We form the matrix

$A:=\begin{pmatrix}-2t_{1}&2\mathbf{a}_{1}^{T}&-1\\\ \vdots&\vdots&\vdots\\\
-2t_{m}&2\mathbf{a}_{m}^{T}&-1\end{pmatrix}\in{\mathbb{R}}^{m\times 5},$ (1.4)

which contains only known quantities. With this, (1.3) can be expressed as a
system of linear equations for the unknown quantities:

$A\cdot\begin{pmatrix}t\\\ \mathbf{x}\\\
\lVert\mathbf{x}\rVert^{2}-t^{2}\end{pmatrix}=\begin{pmatrix}\lVert\mathbf{a}_{1}\rVert^{2}-t_{1}^{2}\\\
\vdots\\\ \lVert\mathbf{a}_{m}\rVert^{2}-t_{m}^{2}\end{pmatrix}.$ (1.5)

Then this has the same solutions as (1.2). Now we make the assumption that $A$
has rank $5$, so (1.5) has a unique solution. (The existence of a solution
follows from the fact that the point $\mathbf{x}$ and time $t$ of emission
exist.) In the case $m=5$ we can then simply invert $A$. If $m>5$, we could
delete all but $5$ linearly independent equations from (1.5), which would give
an algebraically equivalent system with invertible matrix. But in the real
world there are inaccurate measurements, so it should be wiser to apply the
Moore-Penrose inverse $(A^{T}A)^{-1}A^{T}\in{\mathbb{R}}^{5\times m}$.
Specifically, if $B\in{\mathbb{R}}^{4\times m}$ is obtained by deleting the
last row from $(A^{T}A)^{-1}A^{T}$, then

$\begin{pmatrix}t\\\
\mathbf{x}\end{pmatrix}=B\cdot\begin{pmatrix}\lVert\mathbf{a}_{1}\rVert^{2}-t_{1}^{2}\\\
\vdots\\\ \lVert\mathbf{a}_{m}\rVert^{2}-t_{m}^{2}\end{pmatrix}.$ (1.6)

So we have obtained a unique solution for the emission time and place. This is
our first solution, which is available only if $A$ has rank $5$. In 1.4 we
will say something about how likely this it.

But first we will consider the case that $A$ has rank $<5$, and derive our
second solution. What we do assume is that the $\mathbf{a}_{i}$ are not
coplanar. This makes sense, since if the $\mathbf{a}_{i}$ all lay in the same
plane, then even with a known emission time $t$ the location $\mathbf{x}$ of
the source could not be distinguished from the point obtained by reflecting
$\mathbf{x}$ at this plane. Our assumption amounts to saying that the matrix

$\widetilde{A}:=\begin{pmatrix}2\mathbf{a}_{1}^{T}&-1\\\ \vdots&\vdots\\\
2\mathbf{a}_{m}^{T}&-1\end{pmatrix}\in{\mathbb{R}}^{m\times 4}$ (1.7)

has rank $4$ (see Assumption A in [?Beck:Pan:2012]), so in particular we need
$m\geq 4$. The Moore-Penrose inverse is
$\widetilde{B}:=(\widetilde{A}^{T}\widetilde{A})^{-1}\widetilde{A}^{T}\in{\mathbb{R}}^{4\times
m}$, so $\widetilde{B}\widetilde{A}=I_{4}$. Now (1.3) can be restated as

$\widetilde{A}\cdot\begin{pmatrix}\mathbf{x}\\\
\lVert\mathbf{x}\rVert^{2}-t^{2}\end{pmatrix}=2t\begin{pmatrix}t_{1}\\\
\vdots\\\
t_{m}\end{pmatrix}+\begin{pmatrix}\lVert\mathbf{a}_{1}\rVert^{2}-t_{1}^{2}\\\
\vdots\\\ \lVert\mathbf{a}_{m}\rVert^{2}-t_{m}^{2}\end{pmatrix},$ (1.8)

and multiplying by $\widetilde{B}$ yields

$\begin{pmatrix}\mathbf{x}\\\
\lVert\mathbf{x}\rVert^{2}-t^{2}\end{pmatrix}=t\cdot\begin{pmatrix}\mathbf{u}\\\
\alpha\end{pmatrix}+\begin{pmatrix}\mathbf{v}\\\ \beta\end{pmatrix},\
\text{where}\ \begin{pmatrix}\mathbf{u}\\\
\alpha\end{pmatrix}:=2\widetilde{B}\begin{pmatrix}t_{1}\\\ \vdots\\\
t_{m}\end{pmatrix}\ \text{and}\ \begin{pmatrix}\mathbf{v}\\\
\beta\end{pmatrix}:=\widetilde{B}\begin{pmatrix}\lVert\mathbf{a}_{1}\rVert^{2}-t_{1}^{2}\\\
\vdots\\\ \lVert\mathbf{a}_{m}\rVert^{2}-t_{m}^{2}\end{pmatrix}.$ (1.9)

Extracting components, we obtain the equivalent equations

$\mathbf{x}=t\mathbf{u}+\mathbf{v}\quad\text{and}\quad\bigl{(}\lVert\mathbf{u}\rVert^{2}-1\bigr{)}t^{2}+\bigl{(}2\mathbf{u}^{T}\mathbf{v}-\alpha\bigr{)}t+\lVert\mathbf{v}\rVert^{2}-\beta=0.$
(1.10)

Observe that $\mathbf{u}$, $\mathbf{v}$, $\alpha$, and $\beta$ are all derived
from known quantities, so (1.10) can be resolved. In the following theorem,
part (a) summarizes our result in the rank-$5$ case, (b) tells us that (1.10)
is actually equivalent to (1.2), and (c) says that the quadratic equation in
(1.10) never degenerates.

######  Theorem 1.1.

In the above situation and with the notation introduced, we have:

  1. (a)

If the matrix $A$ from (1.4) has rank $5$, then (1.2) has a unique solution
$\left(\begin{smallmatrix}t\\\ \mathbf{x}\end{smallmatrix}\right)$ given by
(1.6).

  2. (b)

Assume that $\operatorname{rank}(A)\leq 4$ and that the $\mathbf{a}_{i}$ are
not coplanar, so in particular $m\geq 4$. Then the equations (1.2) are
satisfied by the same $\left(\begin{smallmatrix}t\\\
\mathbf{x}\end{smallmatrix}\right)$ as (1.10).

  3. (c)

Moreover, the coefficients of $t^{2}$ and $t$ in the quadratic equation in
(1.10) are not both zero, so (1.10) has one or, more likely, two solutions.

###### Proof.

Part (a) has already been shown, so we turn our attention to (b). By
renumbering we may assume that $\mathbf{a}_{1},\ldots,\mathbf{a}_{4}$ are not
coplanar. Then $\widetilde{A}_{4}$, the upper $4\times 4$-part of
$\widetilde{A}$, is invertible, and $A_{4}$, the upper $4\times 5$-part of
$A$, has rank $4$. So we have a matrix $C\in{\mathbb{R}}^{m\times 4}$, with
upper $4\times 4$-part the identity matrix, such that $A=C\cdot A_{4}$ and
$\widetilde{A}=C\cdot\widetilde{A}_{4}$. Extracting the first column from the
first equation gives

$\begin{pmatrix}t_{1}\\\ \vdots\\\
t_{m}\end{pmatrix}=C\cdot\begin{pmatrix}t_{1}\\\ \vdots\\\
t_{4}\end{pmatrix}.$ (1.11)

Since (1.5) has a solution $\left(\begin{smallmatrix}t\\\ \mathbf{x}\\\
\lVert\mathbf{x}\rVert^{2}-t^{2}\end{smallmatrix}\right)$, we obtain

$\begin{pmatrix}\lVert\mathbf{a}_{1}\rVert^{2}-t_{1}^{2}\\\ \vdots\\\
\lVert\mathbf{a}_{m}\rVert^{2}-t_{m}^{2}\end{pmatrix}=A\begin{pmatrix}t\\\
\mathbf{x}\\\
\lVert\mathbf{x}\rVert^{2}-t^{2}\end{pmatrix}=CA_{4}\begin{pmatrix}t\\\
\mathbf{x}\\\
\lVert\mathbf{x}\rVert^{2}-t^{2}\end{pmatrix}=C\cdot\begin{pmatrix}\lVert\mathbf{a}_{1}\rVert^{2}-t_{1}^{2}\\\
\vdots\\\ \lVert\mathbf{a}_{4}\rVert^{2}-t_{4}^{2}\end{pmatrix}.$ (1.12)

We have

$\widetilde{B}=(\widetilde{A}^{T}A)^{-1}\widetilde{A}^{T}=\bigl{(}\widetilde{A}_{4}^{T}C^{T}C\widetilde{A}_{4}\bigr{)}^{-1}\widetilde{A}_{4}^{T}C^{T}=\widetilde{A}_{4}^{-1}(C^{T}C)^{-1}C^{T},$

so

$\widetilde{A}_{4}\widetilde{B}C=I_{4}\quad\text{and}\quad\widetilde{A}\widetilde{B}C=C\widetilde{A}_{4}\widetilde{B}C=C.$
(1.13)

To prove (b), let $\left(\begin{smallmatrix}t\\\
\mathbf{x}\end{smallmatrix}\right)$ satisfy (1.10). Then it satisfies (1.9),
so

$\begin{pmatrix}\mathbf{x}\\\
\lVert\mathbf{x}\rVert^{2}-t^{2}\end{pmatrix}=2t\widetilde{B}\begin{pmatrix}t_{1}\\\
\vdots\\\
t_{m}\end{pmatrix}+\widetilde{B}\begin{pmatrix}\lVert\mathbf{a}_{1}\rVert^{2}-t_{1}^{2}\\\
\vdots\\\
\lVert\mathbf{a}_{m}\rVert^{2}-t_{m}^{2}\end{pmatrix}=2t\widetilde{B}C\begin{pmatrix}t_{1}\\\
\vdots\\\
t_{4}\end{pmatrix}+\widetilde{B}C\begin{pmatrix}\lVert\mathbf{a}_{1}\rVert^{2}-t_{1}^{2}\\\
\vdots\\\ \lVert\mathbf{a}_{4}\rVert^{2}-t_{4}^{2}\end{pmatrix},$

where we used (1.11) and (1.12). With (1.13) this implies

$\widetilde{A}\begin{pmatrix}\mathbf{x}\\\
\lVert\mathbf{x}\rVert^{2}-t^{2}\end{pmatrix}=2tC\begin{pmatrix}t_{1}\\\
\vdots\\\
t_{4}\end{pmatrix}+C\begin{pmatrix}\lVert\mathbf{a}_{1}\rVert^{2}-t_{1}^{2}\\\
\vdots\\\
\lVert\mathbf{a}_{4}\rVert^{2}-t_{4}^{2}\end{pmatrix}=2t\begin{pmatrix}t_{1}\\\
\vdots\\\
t_{m}\end{pmatrix}+\begin{pmatrix}\lVert\mathbf{a}_{1}\rVert^{2}-t_{1}^{2}\\\
\vdots\\\ \lVert\mathbf{a}_{m}\rVert^{2}-t_{m}^{2}\end{pmatrix}.$

Thus (1.8) holds, which is equivalent to (1.2). So every solution of (1.10)
satisfies (1.2), and the converse is true since (1.10) was derived from (1.2).

To prove part (c), assume that the coefficients of $t^{2}$ and $t$ in the
quadratic equation in (1.10) are both zero. Since (1.2) has a solution, this
means that also the constant coefficient vanishes, so there is a solution of
(1.10), and therefore of (1.2), for every $t$. In particular, we have
$\mathbf{x}_{0}\in{\mathbb{R}}^{3}$ such that $\left(\begin{smallmatrix}0\\\
\mathbf{x}_{0}\end{smallmatrix}\right)$ satisfies (1.2). Substituting each
$\mathbf{a}_{i}$ by $\mathbf{a}_{i}-\mathbf{x}_{0}$ preserves the
noncoplanarity of $\mathbf{a}_{1},\ldots,\mathbf{a}_{4}$ and yields new
systems (1.2) and (1.10) of equations. For each solution
$\left(\begin{smallmatrix}t\\\ \mathbf{x}\end{smallmatrix}\right)$ of the old
system (1.2), the new one now has the solution $\left(\begin{smallmatrix}t\\\
\mathbf{x}-\mathbf{x}_{0}\end{smallmatrix}\right)$, so the new quadratic
equation in (1.10) still has infinitely many solutions. Moreover,
$\left(\begin{smallmatrix}0\\\ \mathbf{0}\end{smallmatrix}\right)$ is a
solution of (1.2), so $\lVert\mathbf{a}_{i}\rVert=|t_{i}|$. Hence (1.9)
implies $\mathbf{v}=\mathbf{0}$ and $\beta=0$. The vanishing of the
coefficient of $t$ in the quadratic equation therefore means $\alpha=0$,
therefore means $\alpha=0$, so

$2\begin{pmatrix}\mathbf{a}_{1}^{T}\mathbf{u}\\\ \vdots\\\
\mathbf{a}_{4}^{T}\mathbf{u}\end{pmatrix}=\begin{pmatrix}2\mathbf{a}_{1}^{T}&-1\\\
\vdots&\vdots\\\
2\mathbf{a}_{4}^{T}&-1\end{pmatrix}\cdot\begin{pmatrix}\mathbf{u}\\\
\alpha\end{pmatrix}=\widetilde{A}_{4}\begin{pmatrix}\mathbf{u}\\\
\alpha\end{pmatrix}\underset{\eqref{1eqUV}}{=}\\\
2\widetilde{A}_{4}\widetilde{B}\begin{pmatrix}t_{1}\\\ \vdots\\\
t_{m}\end{pmatrix}\underset{\eqref{1eqT}}{=}2\widetilde{A}_{4}\widetilde{B}C\begin{pmatrix}t_{1}\\\
\vdots\\\
t_{4}\end{pmatrix}\underset{\eqref{1eqMP}}{=}2\begin{pmatrix}t_{1}\\\
\vdots\\\
t_{4}\end{pmatrix}=2\begin{pmatrix}\varepsilon_{1}\lVert\mathbf{a}_{1}\rVert\\\
\vdots\\\ \varepsilon_{4}\lVert\mathbf{a}_{4}\rVert\end{pmatrix}$

for some $\varepsilon_{i}\in\\{\pm 1\\}$. The vanishing of the
$t^{2}$-coefficient in (1.10) means that $\lVert u\rVert=1$, so

$|\mathbf{a}_{i}^{T}\mathbf{u}|=\lVert\mathbf{a}_{i}\rVert=\lVert\mathbf{a}_{i}\rVert\cdot\lVert
u\rVert\quad(i=1,\ldots,4).$

Thus the Cauchy–Schwarz inequality is actually an equality, implying that
$\mathbf{a}_{i}$ and $\mathbf{u}$ are linearly dependent. This shows that
$\mathbf{a}_{1},\ldots,\mathbf{a}_{4}$ are collinear, contradicting the
hypothesis that they are not coplanar. ∎

######  Remark 1.2.

Everything in this section carries over directly to $n$-dimensional space,
with $n>1$. Just replace every instance of $3$, $4$, and $5$ by $n$, $n+1$,
and $n+2$, and replace “coplanar” by “contained in a common
$(n-1)$-dimensional affine subspace.”

What happens for $n=1$? Everything works, except for the very last sentence in
the proof of Theorem 1.1: In dimension $1$, collinearity does not contradict
being noncopunctual (i.e., not being the same point). But this makes
everything break down. In fact, the mathematics bears out what has always been
clear about the one-dimensional case: If the source lies on the same side of
every sensor, then there is no way to find out its position from the time
differences of signal arrivals; and indeed in this case the coefficients of
the quadratic equation in (1.10) are all $0$, and the matrix
$A\in{\mathbb{R}}^{m\times 3}$ has rank $2$.

Even though the one-dimensional case is not interesting in itself, it shows,
as do (5) and (6) in Example 1.3 below, that the effort of proving Theorem
1.1(c) was not irrelevant: this is not a truism. $\triangleleft$

Up to now, we have worked with the equations (1.2), and seen that they may
have two solutions. But according to (1.1) (which expresses that signals
arrive after having been sent) we have $t_{i}\geq t$. If one of the solutions
does not satisfy this, it is spurious and can be discarded. But if both do,
the given data do not uniquely determine $\mathbf{x}$ and $t$. The first two
of the following examples show that this can actually happen. The third
example has a spurious solution, and the others exemplify some special cases.

######  Example 1.3.

In the following, we chose the coordinates in a way to keep all numbers
rational, so the examples are quite Pythagorean-triple-prone.

  1. (1)

Of the five sensor positions

$\mathbf{a}_{1}=\begin{pmatrix}3\\\ 4\\\ 0\end{pmatrix},\
\mathbf{a}_{2}=\begin{pmatrix}-2\\\ -2\\\ 1\end{pmatrix},\
\mathbf{a}_{3}=\begin{pmatrix}-1\\\ 0\\\ 0\end{pmatrix},\
\mathbf{a}_{4}=\frac{2}{21}\begin{pmatrix}0\\\ -24\\\ 7\end{pmatrix},\
\mathbf{a}_{5}=\frac{1}{21}\begin{pmatrix}0\\\ 76\\\ 0\end{pmatrix},$

no four are coplanar. If we assume that the source is at the origin
$\mathbf{x}=\mathbf{0}=(0,0,0)^{T}$ and emits a signal at time $t=0$, then the
$i$-th sensor will receive this at time $t_{i}=\lVert\mathbf{a}_{i}\rVert$. We
have $(t_{1},t_{2},t_{3},t_{4},t_{5})=(5,3,1,50/21,76/21)$. We chose the
$\mathbf{a}_{i}$ in such a way that the affine relation
$2\mathbf{a}_{1}+2\mathbf{a}_{2}+2\mathbf{a}_{3}-3\mathbf{a}_{4}-3\mathbf{a}_{5}=\mathbf{0}$
is also satisfied by their lengths $\lVert\mathbf{a}_{i}\rVert$. So the matrix
$A\in{\mathbb{R}}^{5\times 5}$ is not invertible. The computations according
to (1.9) yield

$\mathbf{u}=\frac{1}{55}\begin{pmatrix}21\\\ 34\\\ 199\end{pmatrix},\
\mathbf{v}=\mathbf{0},\ \alpha=-\frac{152}{55},\ \text{and}\ \beta=0,$

so (1.10) becomes

$\mathbf{x}=t\mathbf{u}\quad\text{and}\quad\frac{38173}{3025}t^{2}+\frac{152}{55}t=0,$

which is solved by

$t=0,\ \mathbf{x}=\mathbf{0},\quad\text{and}\quad
t^{\prime}=\frac{-8360}{38173},\
\mathbf{x}^{\prime}=\frac{-152}{38173}\begin{pmatrix}21\\\ 34\\\
199\end{pmatrix}.$

Now it can easily be verified directly that (1.1) holds with $t$ and
$\mathbf{x}$ replaced by $t^{\prime}$ and $\mathbf{x}^{\prime}$. In other
words, had the signal been emitted from the position $\mathbf{x}^{\prime}$ at
time $t^{\prime}$ rather than from $\mathbf{x}=\mathbf{0}$ at $t=0$, it would
have arrived at the exact same times $t_{i}$ at the sensors. So even if the
system (1.1) is overdetermined ($5$ equations for $4$ unknowns) and the sensor
positions are not chosen in an obviously clumsy way, it may still be
impossible to uniquely determine the source position.

  2. (2)

A simpler example of the same type can be constructed in $2$ dimensions. Take
the sensor positions

$\mathbf{a}_{1}=\begin{pmatrix}9\\\ 12\end{pmatrix},\
\mathbf{a}_{2}=\begin{pmatrix}9\\\ -12\end{pmatrix},\
\mathbf{a}_{3}=\begin{pmatrix}10\\\ -24\end{pmatrix},\
\mathbf{a}_{4}=\begin{pmatrix}10\\\ 24\end{pmatrix}.$

The source is again located at $\mathbf{x}=\mathbf{0}$ and emits its signal at
time $t=0$. It is then received at times
$(t_{1},t_{2},t_{3},t_{4})=(15,15,26,26)$. The same computation as above (but
with smaller numbers) shows that emission time and place

$t^{\prime}=\frac{7}{5}\quad\text{and}\quad\mathbf{x}^{\prime}=\frac{77}{5}\begin{pmatrix}1\\\
0\end{pmatrix}$

would have led to the exact same arrival times $t_{i}$.

  3. (3)

An example with a spurious solution (again in dimension $2$) is given by

$\mathbf{a}_{1}=\begin{pmatrix}4\\\ 0\end{pmatrix},\
\mathbf{a}_{2}=\begin{pmatrix}-3\\\ 4\end{pmatrix},\
\mathbf{a}_{3}=\begin{pmatrix}-3\\\ -4\end{pmatrix},\
\mathbf{x}=\begin{pmatrix}0\\\ 0\end{pmatrix},\ \text{and}\ t=0.$

Reception times are $(t_{1},t_{2},t_{3})=(4,5,5)$. Here another solution of
(1.2) is $t^{\prime}=28/3$ and
$\mathbf{x}^{\prime}=-\frac{4}{3}\left(\begin{smallmatrix}1\\\
0\end{smallmatrix}\right)$. This is spurious, since $t_{i}-t^{\prime}<0$. In
fact, if the signal had been sent from $\mathbf{x}^{\prime}$ at time
$t^{\prime}$, it would have arrived at times
$(\frac{44}{3},\frac{41}{3},\frac{41}{3})$. The absolute differences between
the arrival times are the same, but the sequence is reversed.

  4. (4)

It can happen that the signal arrives at the same time $t_{1}=\cdots=t_{m}$ at
all sensors, so the source has the same distance from all of them. Intuition
tells us that there can only be one such point, and again the mathematics
bears this out. In fact, the equation $\widetilde{B}\widetilde{A}=I_{n+1}$
before (1.8) and the definition of $\widetilde{A}$ imply that all rows from
$\widetilde{B}$, except for the last one, have coefficient sum $0$. Therefore
$\mathbf{u}=0$, so indeed there is only one solution for $\mathbf{x}$. Of the
solutions for $t$, one is spurious.

  5. (5)

Assume that $\mathbf{x}$ is collinear with two of the $\mathbf{a}_{i}$, say
$\mathbf{a}_{1}$ and $\mathbf{a}_{2}$, but does not lie between them. Choosing
the coordinate system suitably, we may assume $\mathbf{x}=\mathbf{0}$ and
$t=0$. Then our assumption means $\mathbf{a}_{1}=\lambda\mathbf{a}_{2}$ with
$1\neq\lambda>0$. Now $\widetilde{B}\widetilde{A}=I_{n+1}$ implies that the
last row of $\widetilde{B}$ is
$(\frac{-1}{\lambda-1},\frac{\lambda}{\lambda-1},0,\ldots,0)$. We have
$t_{1}=\lVert\mathbf{a}_{1}\rVert=\lambda\lVert\mathbf{a}_{2}\rVert=\lambda
t_{2}$, so $\alpha=0$ by (1.9). Since $\mathbf{v}=\mathbf{0}$ and $\beta=0$,
the quadratic equation becomes $t^{2}=0$.

  6. (6)

Here is an example where the coefficient of $t^{2}$ becomes $0$:

$\mathbf{a}_{1}=\begin{pmatrix}1\\\ 0\end{pmatrix},\
\mathbf{a}_{2}=\begin{pmatrix}-1\\\ 0\end{pmatrix},\
\mathbf{a}_{3}=\begin{pmatrix}3\\\ 4\end{pmatrix},\
\mathbf{x}=\begin{pmatrix}0\\\ 0\end{pmatrix},\ \text{and}\ t=0.$

The computation shows $\mathbf{u}=\left(\begin{smallmatrix}0\\\
1\end{smallmatrix}\right)$, $\mathbf{v}=\left(\begin{smallmatrix}0\\\
0\end{smallmatrix}\right)$, $\alpha=-2$, and $\beta=0$. So here the quadratic
equation degenerates to $2t=0$. This can be interpreted as follows: As the
positions $\mathbf{a}_{i}$ approach the values given above, the alternative
solution $\left(\begin{smallmatrix}t^{\prime}\\\
\mathbf{x}^{\prime}\end{smallmatrix}\right)$ of (1.10), apart from the
solution $\left(\begin{smallmatrix}t\\\
\mathbf{x}\end{smallmatrix}\right)=\left(\begin{smallmatrix}0\\\
\mathbf{0}\end{smallmatrix}\right)$, tends to infinity. For example, taking
$\mathbf{a}_{3}=\left(\begin{smallmatrix}3\\\ 3.99\end{smallmatrix}\right)$
and leaving $\mathbf{a}_{1}$ and $\mathbf{a}_{2}$ unchanged leads to
$t^{\prime}\approx-1991$ and
$\mathbf{x}^{\prime}\approx\left(\begin{smallmatrix}0\\\
-1992\end{smallmatrix}\right)$. $\triangleleft$

Having seen from Example 1.3(1) that even with $5$ sensors such that no $4$ of
them are coplanar, it may happen that the matrix $A$, defined in (1.4) is not
invertible, we wonder how often this happens. The following result says that
under mild hypotheses, the answer is “very rarely”, i.e., almost certainly the
formula (1.6) can be applied for finding $\mathbf{x}$ and $t$. Notice that $A$
is formed with times $t_{i}$ given by (1.1), so, because of the last column of
$A$, its rank only depends on the positions $\mathbf{a}_{i}$ and $\mathbf{x}$.

######  Proposition 1.4.

Assume we have $m=5$ sensors such that

$\det\begin{pmatrix}\varepsilon_{1}\cdot\lVert\mathbf{a}_{1}\rVert&\mathbf{a}_{1}^{T}&1\\\
\vdots&\vdots&\vdots\\\
\varepsilon_{5}\cdot\lVert\mathbf{a}_{5}\rVert&\mathbf{a}_{5}^{T}&1\end{pmatrix}\neq
0\quad\text{for all}\quad\varepsilon_{1},\ldots,\varepsilon_{5}=\pm 1.$ (1.14)

Then the set of all $\mathbf{x}\in{\mathbb{R}}^{3}$ such that the matrix $A$
has rank $<5$ is contained in a $2$-dimensional subvariety of
${\mathbb{R}}^{3}$.

###### Proof.

The function

$f(\mathbf{x},\mathbf{a}_{1},\ldots,\mathbf{a}_{5}):=\prod_{\varepsilon_{1},\ldots,\varepsilon_{5}=\pm
1}\det\begin{pmatrix}\varepsilon_{1}\cdot\lVert\mathbf{a}_{1}-\mathbf{x}\rVert&\mathbf{a}_{1}^{T}&1\\\
\vdots&\vdots&\vdots\\\
\varepsilon_{5}\cdot\lVert\mathbf{a}_{5}-\mathbf{x}\rVert&\mathbf{a}_{5}^{T}&1\end{pmatrix}$

is a polynomial in the coefficients of $\mathbf{x}$ and the $\mathbf{a}_{i}$.
If $A$ has rank $<5$, then
$f(\mathbf{x},\mathbf{a}_{1},\ldots,\mathbf{a}_{5})=0$. Our hypothesis means
that $f(\mathbf{0},\mathbf{a}_{1},\ldots,\mathbf{a}_{5})\neq 0$, so the
assertion follows. ∎

The hypothesis (1.14) in 1.4 looks a bit messy and lacks geometric content,
but is readily verifiable. It would be desirable to have some more geometric
conditions under which the assertion of 1.4 holds. The following conjecture
would be the best possible result, since its converse is clearly true.

######  Conjecture 1.5.

The assertion of 1.4 holds under the milder hypothesis that the
$\mathbf{a}_{i}$ are not coplanar and pairwise distinct.

We managed to prove the conjecture in the $2$-dimensional case by considering
the polynomial $f$ used in the above proof as a polynomial in the coordinates
of $\mathbf{x}$ as main variables, and forming the ideal generated by the
coefficients. A computation in the computer algebra system MAGMA [?magma] then
shows that the equations that express that the $\mathbf{a}_{i}$ do not satisfy
the hypothesis of the conjecture all lie in the radical of this ideal.
However, in the $3$-dimensional case our computations ran into an impasse. The
$1$-dimensional case of the conjecture is false.

######  Remark 1.6.

The first, and most cited, algebraic solution of the GPS-problem appears to
have been given by ?Bancroft:1985 [?Bancroft:1985] in ?Bancroft:1985. Let us
point out some differences between his approach and ours.

  * •

Bancroft reaches a quadratic equation even if there are more than $4$ sensors,
and does not offer an explicit formula such as (1.6).

  * •

For reaching the equations (1.10) in the case that $A$ has rank $<5$, we
assume that $\widetilde{A}$, defined in (1.7), has rank $4$ or, equivalently,
that the $\mathbf{a}_{i}$ are not coplanar. On the other hand, Bancroft
assumes that the matrix

$\begin{pmatrix}\mathbf{a}_{1}^{T}&t_{1}\\\ \vdots&\vdots\\\
\mathbf{a}_{m}^{T}&t_{m}\end{pmatrix}\in{\mathbb{R}}^{m\times 4}$

has rank $4$. But whether this is the case does not only depend on the
positions $\mathbf{a}_{i}$ and $\mathbf{x}$, but also on $t$, which in
Bancroft’s situation is the bias between the clocks. For example, if there are
$m=4$ noncoplanar sensors, there is always a value for $t$ such that the above
matrix has determinant $0$.

  * •

There is no proof that Bancroft’s quadratic equation does not degenerate.
$\triangleleft$

##  2\. Relations and matching

We consider the same situation as in the previous section: A source at
position $\mathbf{x}\in{\mathbb{R}}^{3}$ emits a signal at time $t$, which is
received by $m$ sensors at positions
$\mathbf{a}_{1},\ldots,\mathbf{a}_{m}\in{\mathbb{R}}^{3}$ and at times
$t_{1},\ldots,t_{m}$. The unit of time is chosen such that the signal
propagates with speed $1$, so $\lVert\mathbf{a}_{i}-\mathbf{x}\rVert=t_{i}-t$.
The $\mathbf{a}_{i}$ and $t_{i}$ are considered as known, and Section 1 was
about how to find $\mathbf{x}$ and $t$ from them. Now if there are $5$ sensors
or more, then, according to the following result, there are algebraic
relations between the known quantities.

######  Theorem 2.1.

In the above situation, write
$d_{i,j}:=\lVert\mathbf{a}_{i}-\mathbf{a}_{j}\rVert$ and
$t_{i,j}:=t_{i}-t_{j}$. Then the matrix

$D=\bigl{(}t_{i,j}^{2}-d_{i,j}^{2}\bigr{)}_{i,j=1,\ldots,m}=\begin{pmatrix}0&t_{1,2}^{2}-d_{1,2}^{2}&\cdots&t_{1,m}^{2}-d_{1,m}^{2}\\\
t_{2,1}^{2}-d_{2,1}^{2}&0&\cdots&t_{2,m}^{2}-d_{2,m}^{2}\\\
\vdots&\vdots&\ddots&\vdots\\\
t_{m,1}^{2}-d_{m,1}^{2}&t_{m,2}^{2}-d_{m,2}^{2}&\cdots&0\end{pmatrix}\in{\mathbb{R}}^{m\times
m}$

has rank $\leq 4$. So if $m\geq 5$ we have the relation $\det(D)=0$, and if
$m>5$ all $5\times 5$-minors are zero.

We have formulated Theorem 2.1 in the $3$-dimensional case for the sake of
simplicity. But it actually holds in any dimension $n$, saying that
$\operatorname{rank}(D)\leq n+1$.

For the proof we will use the Cayley-Menger matrix. The following result about
its rank is well known (see ?Cayley:1841 [?Cayley:1841]) in the case of
Euclidean spaces, but we need a more general version given by the following
proposition. We will not need the exact value of the rank, but include it for
the sake of completeness.

######  Proposition 2.2.

Let $\mathbf{v}_{0},\ldots,\mathbf{v}_{m}\in V$ be vectors in a Euclidean
space or, more generally, in a vector space over a field of characteristic
$\neq 2$ equipped with a quadratic form $q$. Set
$\delta_{i,j}:=q(\mathbf{v}_{i}-\mathbf{v}_{j})$, which in the special case of
a Euclidean space is the squared distance between $\mathbf{v}_{i}$ and
$\mathbf{v}_{j}$. Then the Cayley-Menger matrix

$C:=\begin{pmatrix}0&1&1&1&\cdots&1\\\
1&0&\delta_{0,1}&\delta_{0,2}&\cdots&\delta_{0,m}\\\
1&\delta_{1,0}&0&\delta_{1,2}&\cdots&\delta_{1,m}\\\
1&\delta_{2,0}&\delta_{2,1}&0&\cdots&\delta_{2,m}\\\
\vdots&\vdots&\vdots&\vdots&\ddots&\vdots\\\
1&\delta_{m,0}&\delta_{m,1}&\delta_{m,2}&\cdots&0\end{pmatrix}\in
K^{(m+2)\times(m+2)}$

has rank $\leq\dim(V)+2$. More precisely, if $r$ is the rank of $q$ restricted
to the subspace $U\subseteq V$ generated by
$\mathbf{v}_{1}-\mathbf{v}_{0},\ldots,\mathbf{v}_{m}-\mathbf{v}_{0}$, then
$\operatorname{rank}(C)=r+2$. Notice that $r\leq\dim(U)$, with equality in the
special case of a Euclidean space. Also notice that $\dim(U)$ is equal to the
dimension of the affine subspace generated by
$\mathbf{v}_{0},\ldots,\mathbf{v}_{m}$.

###### Proof.

We only need to show $\operatorname{rank}(C)=r+2$ in the more general case of
a quadratic space over a field $K$. Replacing each $\mathbf{v}_{i}$ by
$\mathbf{v}_{i}-\mathbf{v}_{0}$ does not change $C$ or $U$, so by doing this
we may assume $\mathbf{v}_{0}=\mathbf{0}$. Then $U$ is generated (as a vector
space) by $\mathbf{v}_{1},\ldots,\mathbf{v}_{m}$, and we may replace $V$ by
$U$. Now $V$ is generated by the $\mathbf{v}_{i}$ and in particular finite-
dimensional. By choosing a basis of $V$ we may then replace $V$ by $K^{n}$.
Then $q$ is given by $q(\mathbf{v})=\mathbf{v}^{T}A\mathbf{v}$ for
$\mathbf{v}\in V=K^{n}$, with $A\in K^{n\times n}$ a symmetric matrix of rank
$r$. $A$ also defines the bilinear form $\langle\cdot,\cdot\rangle$ belonging
to $q$. After these reductions, the main part of the proof rests on the
following matrix computation.

With

$E:=\left(\begin{array}[]{c|c}\begin{array}[]{cc}0&1\\\ 1&0\end{array}&0\\\
\hline\cr\vphantom{3^{3^{3^{3}}}}0&-2A\end{array}\right)\in
K^{(n+2)\times(n+2)},\
F:=\left(\begin{array}[]{c|c}I_{2}&\begin{array}[]{ccc}\delta_{0,1}&\cdots&\delta_{0,m}\\\
1&\cdots&1\end{array}\\\
\hline\cr\vphantom{3^{3^{3^{3}}}}0&\begin{array}[]{ccc}\mathbf{v}_{1}&\cdots&\mathbf{v}_{m}\end{array}\end{array}\right)\in
K^{(n+2)\times(m+2)}$

we have

$F^{T}EF=\left(\begin{array}[]{c|c}I_{2}&\vphantom{3_{3_{3_{3}}}}0\\\
\hline\cr\begin{array}[]{cc}\delta_{0,1}&\vphantom{3^{3^{3^{3}}}}1\\\
\vdots&\vdots\\\
\delta_{0,m}&1\end{array}&\begin{array}[]{c}\mathbf{v}_{1}^{T}\\\ \vdots\\\
\mathbf{v}_{m}^{T}\end{array}\end{array}\right)\left(\begin{array}[]{c|c}\begin{array}[]{cc}0&1\\\
1&0\end{array}&\begin{array}[]{ccc}1&\cdots&\hphantom{A}1\\\
\hphantom{a}\delta_{0,1}&\cdots&\hphantom{-2}\delta_{0,m}\end{array}\\\
\hline\cr\vphantom{3^{3^{3^{3}}}}0&\begin{array}[]{ccc}-2A\mathbf{v}_{1}&\cdots&-2A\mathbf{v}_{m}\end{array}\end{array}\right)\\\
=\left(\begin{array}[]{c|c}\begin{array}[]{cc}0&\hphantom{2}1\hphantom{2}\\\
1&\hphantom{2}0\hphantom{2}\end{array}&\begin{array}[]{ccc}1&\cdots&1\\\
\delta_{0,1}&\cdots&\delta_{0,m}\end{array}\\\
\hline\cr\begin{array}[]{cc}\vphantom{3^{3^{3^{3}}}}1&\delta_{1,0}\\\
\vdots&\vdots\\\ 1&\delta_{m,0}\end{array}&B\end{array}\right),$

where the $(i,j)$-th entry of the matrix $B$ is

$\delta_{0,i}+\delta_{0,j}-2\mathbf{v}_{i}^{T}A\mathbf{v}_{j}=\langle\mathbf{v}_{i},\mathbf{v}_{i}\rangle+\langle\mathbf{v}_{j},\mathbf{v}_{j}\rangle-2\langle\mathbf{v}_{i},\mathbf{v}_{j}\rangle=\langle\mathbf{v}_{i}-\mathbf{v}_{j},\mathbf{v}_{i}-\mathbf{v}_{j}\rangle=\delta_{i,j}.$

So $F^{T}EF=C$. Since the $\mathbf{v}_{i}$ span $V=K^{n}$, $F$ has rank $n+2$,
so the linear map $K^{m+2}\to K^{n+2}$ given by $F$ is surjective. Likewise,
the map given by $F^{T}$ is injective, and the image of the map given by $E$
has dimension equal to $\operatorname{rank}(A)+2=r+2$. It follows that the map
given by $C$ has an image of dimension $r+2$, which is our claim. ∎

###### Proof of Theorem 2.1.

We apply 2.2 with $V={\mathbb{R}}^{4}$ being Minkowski space, so for
$\mathbf{v}=\left(\begin{smallmatrix}t\\\
\mathbf{u}\end{smallmatrix}\right)\in V$ with $t\in{\mathbb{R}}$,
$\mathbf{u}\in{\mathbb{R}}^{3}$ the quadratic form is
$q(\mathbf{v})=t^{2}-\lVert\mathbf{u}\rVert^{2}$, with $\lVert\cdot\rVert$ the
usual Euclidean norm. For $i=1,\ldots,m$, set
$\mathbf{v}_{i}:=\left(\begin{smallmatrix}t_{i}\\\
\mathbf{a}_{i}\end{smallmatrix}\right)$, and also set
$\mathbf{v}_{0}:=\left(\begin{smallmatrix}t\\\
\mathbf{x}\end{smallmatrix}\right)$ (the time and place of emission). Then the
equation $\lVert\mathbf{a}_{i}-\mathbf{x}\rVert=t_{i}-t$ implies
$\delta_{0,i}=q(\mathbf{v}_{0}-\mathbf{v}_{i})=0$ for all $i$. So the Cayley-
Menger matrix becomes

$C=\left(\begin{array}[]{c|c}\begin{array}[]{cc}0&1\\\
1&0\end{array}&\begin{array}[]{ccc}1&\cdots&1\\\ 0&\cdots&0\end{array}\\\
\hline\cr\begin{array}[]{cc}\vphantom{3^{3^{3}}}1&0\\\ \vdots&\vdots\\\
1&0\end{array}&D\end{array}\right),$

with $D$ as defined in the theorem. Since $\operatorname{rank}(C)\leq 6$ by
2.2, we get $\operatorname{rank}(D)\leq 4$. ∎

By Theorem 2.1 there is a relation between the reception times of a signal
coming from a single emission event. If there are multiple emission events
(from different source locations and/or at different times), this relation can
be used to match those reception times that come from the same event. Matching
reception times can then be used to determine the time and place of emission,
making use of the methods from Section 1. Algorithm 2.3 makes this idea
precise.

Algorithm 2.3 Detect source positions and emission times of multiple emission
events

Input:

For $i=1,\ldots,5$, a set $\mathcal{T}_{i}$ containing the points in time when
the $i$th sensor received a signal. The units of time and distance should be
chosen such that signals travel with speed $1$. The positions
$\mathbf{a}_{1},\ldots,\mathbf{a}_{5}\in{\mathbb{R}}^{3}$ of the sensors need
to be known.

1: Let $d_{i,j}:=\lVert\mathbf{a}_{i}-\mathbf{a}_{j}\rVert$ be the distances
between the sensors. Set $\mathcal{E}:=\emptyset$. The detected emission
events will be collected as pairs $(t,\mathbf{x})$ in the set $\mathcal{E}$.

2: for
$(t_{1},t_{2},t_{3},t_{4},t_{5})\in\mathcal{T}_{1}\times\mathcal{T}_{2}\times\mathcal{T}_{3}\times\mathcal{T}_{4}\times\mathcal{T}_{5}$
do

3: Set up the matrix
$D=\bigl{(}(t_{i}-t_{j})^{2}-d_{i,j}^{2}\bigr{)}_{i,j=1,\ldots,5}$.

4: if $\det(D)=0$ then

5: Use Theorem 1.1, with the current $t_{1},\ldots,t_{5}$ and
$\mathbf{a}_{1},\ldots,\mathbf{a}_{5}$ as input, to compute the emission time
$t$ and the source position $\mathbf{x}$. Include $(t,\mathbf{x})$ in the set
$\mathcal{E}$. In the unlikely event that Theorem 1.1 yields two solutions and
neither can be discarded (as spurious or from other context), include both.

6: end if

7: end for

8: Output the set $\mathcal{E}$ of detected emission events.

######  Remark 2.4.

In fact, Algorithm 2.3 does more than just matching reception times of
signals. It also discards spurious signals registered by sensors. By this we
mean erroneous registrations of signals, or registrations of signals that
originate very near to a sensor and are irrelevant since they cannot be
perceived by other sensors. Indeed, it is almost impossible for such a
spurious signal to satisfy the relation $\det(D)=0$ together with other
”legitimate” reception times. $\triangleleft$

In a later paper we will study the behavior (and modifications) of the
algorithm in situations where the input data is inexact because of measurement
errors.

##  3\. Wall detection by echoes

Theorem 2.1 guarantees that the relation $\det(D)=0$ always holds if the
signals received at times $t_{1},\ldots,t_{5}$ come from the same emission
event. So all events for which a signal is received by every sensor will be
detected. But it is possible that the determinant becomes $0$ even if the
signals do not come from the same event. For instance, if $t_{1},\ldots,t_{4}$
do come from the same event, a different source may send a signal at such an
unlucky time that it is received by the fifth sensor at a time $t_{5}$ that
happens to make the equation $\det(D)=0$ come true. This can happen no matter
where sources and sensors are positioned. In fact, with the signal reception
times as the only available information, any algorithm would be tricked into
making a false match if some signal arrives at an unlucky time.

We will now restrict our attention to the case in which all emission events
share the same emission time (which is still unknown to the sensors). This
happens when the emission events are in fact echoes of a single sound emission
bouncing off from various walls (i.e., flat surfaces). In fact, in the ray
acoustics approximation the received echoes virtually come from the so-called
mirror points, i.e., the points obtained by reflecting the original source
position at the walls. If an echo from a wall is received by all five sensors
(or microphones in the acoustic case), then Algorithm 2.3 computes the mirror
point, and it is easy to find the wall position from this. With the
restriction to simultaneous (virtual) emission events, it becomes less likely
that Algorithm 2.3 produces a mismatch ($\det(D)=0$ even though the signals
come from different mirror points) and so erroneously detects a wall which is
not really there (often called a ghost wall).

So we can be hopeful that, in contrast to the situation with different
emission times, the choice of the sensor positions may preclude ghost walls.
To give this a name, we say that the sensors are in a good position if
Algorithm 2.3 produces no ghost walls. Whether this is true clearly also
depends on the coordinates of the mirror points, but not on the time of the
sound emission, since only differences of reception times go into the matrix
$D$. Even more hopeful, we say that almost all positions are good if the bad
positions are contained in a lower-dimensional subvariety of the configuration
space $({\mathbb{R}}^{3})^{5}$ of all possible microphone positions.
Intuitively “almost all” can be thought of as “with probability one.”

######  Theorem 3.1.

Consider a given room, by which we understand an arrangement of walls, which
may include ceilings, floors, and sloping walls. Assume there is a loudspeaker
at a given position in the room. Now five microphones are positioned in the
room. Then almost all loudspeaker positions are good, meaning that from a
single sound emitted by the loudspeaker, Algorithm 2.3 detects all walls from
which an echo is received by every microphone, but it detects no ghost walls.

###### Proof.

We are given a finite set $\mathcal{S}\subset{\mathbb{R}}^{3}$ of mirror
points, obtained by reflecting the given position of the loudspeaker at the
various walls. From now on we can forget about the walls and the loudspeaker,
since the signals will be received by the microphones as if they were all
simultaneously emitted from the mirror points. The microphone positions can be
represented as the columns of a matrix
$M=(\mathbf{a}_{1},\ldots,\mathbf{a}_{5})\in{\mathbb{R}}^{3\times 5}$. We can
thus speak of good or bad matrices $M$. We will also say that $M$ is very good
if the following is true: For any five points
$\mathbf{s}_{1},\ldots,\mathbf{s}_{5}\in\mathcal{S}$, the relation

$\det\Bigl{(}\bigl{(}\lVert\mathbf{s}_{i}-\mathbf{a}_{i}\rVert-\lVert\mathbf{s}_{j}-\mathbf{a}_{j}\rVert\bigl{)}^{2}-\lVert\mathbf{a}_{i}-\mathbf{a}_{j}\rVert^{2}\Bigr{)}_{i,j=1,\ldots,5}=0$
(3.1)

only holds if
$\mathbf{s}_{1}=\mathbf{s}_{2}=\mathbf{s}_{3}=\mathbf{s}_{4}=\mathbf{s}_{5}$.
A very good position is a good one, since for a quintuple
$(t_{1},\ldots,t_{5})\in\mathcal{T}_{1}\times\cdots\times\mathcal{T}_{5}$
considered in Algorithm 2.3 there are points $\mathbf{s}_{i}\in\mathcal{S}$
such that $t_{i}-t=\lVert\mathbf{s}_{i}-\mathbf{a}_{i}\rVert$ (with $t$ the
time of simultaneous emission). So the determinant considered in the algorithm
is just the one in (3.1). If this is zero, the hypothesis of “very good”
implies that the $\mathbf{s}_{i}$ are all the same point
$\mathbf{s}\in\mathcal{S}$. So indeed if the algorithm includes a point into
the set $\mathcal{E}$, this will be an actual mirror point, meaning that $M$
is good.

The main difficulty of the proof arises from the fact that the determinant in
(3.1) is not a polynomial in the coordinates of the $\mathbf{a}_{i}$ and the
$\mathbf{s}_{i}$, because the norms involve square roots. To circumvent this
problem, we form the product

$f(\mathbf{a}_{1},\ldots,\mathbf{a}_{5},\mathbf{s}_{1},\ldots,\mathbf{s}_{5}):=\\!\\!\\!\\!\prod_{\begin{subarray}{c}\varepsilon_{1},\ldots,\varepsilon_{4}=\pm
1\\\
\epsilon_{5}=1\end{subarray}}\\!\\!\\!\\!\det\Bigl{(}\bigl{(}\varepsilon_{i}\lVert\mathbf{s}_{i}-\mathbf{a}_{i}\rVert-\varepsilon_{j}\lVert\mathbf{s}_{j}-\mathbf{a}_{j}\rVert\bigl{)}^{2}-\lVert\mathbf{a}_{i}-\mathbf{a}_{j}\rVert^{2}\Bigr{)}_{i,j=1,\ldots,5},$

which is easily seen to be a polynomial in the coordinates of its arguments.
Let us say that $M$ is excellent if for any
$\mathbf{s}_{1},\ldots,\mathbf{s}_{5}\in\mathcal{S}$ the relation
$f(\mathbf{a}_{1},\ldots,\mathbf{a}_{5},\mathbf{s}_{1},\ldots,\mathbf{s}_{5})=0$
implies that the $\mathbf{s}_{i}$ are all equal. So “excellent” implies “very
good” and “good”. Now the set

$\mathcal{U}_{\mathbf{s}_{1},\ldots,\mathbf{s}_{5}}:=\bigl{\\{}M=(\mathbf{a}_{1},\ldots,\mathbf{a}_{5})\in{\mathbb{R}}^{3\times
5}\mid
f(\mathbf{a}_{1},\ldots,\mathbf{a}_{5},\mathbf{s}_{1},\ldots,\mathbf{s}_{5})\neq
0\bigr{\\}}$

is Zariski open in ${\mathbb{R}}^{3\times 5}$, and

$\mathcal{U}:=\bigcap_{\begin{subarray}{c}\mathbf{s}_{1},\ldots,\mathbf{s}_{5}\in
S\ \text{such that}\\\ \text{not all}\ \mathbf{s}_{i}\ \text{are
equal}\end{subarray}}\mathcal{U}_{\mathbf{s}_{1},\ldots,\mathbf{s}_{5}}$

is the set of excellent matrices. So if we can show that
$\mathcal{U}_{\mathbf{s}_{1},\ldots,\mathbf{s}_{5}}\neq\emptyset$ for all
$\mathbf{s}_{1},\ldots,\mathbf{s}_{5}\in\mathcal{S}$ that are not all equal,
by the Zariski openness the theorem follows. Since we have no control over the
given set $\mathcal{S}$, we need to show the nonemptiness for any five points
$\mathbf{s}_{i}\in{\mathbb{R}}^{3}$ that are not all equal, and this also
suffices. Equivalently, we need to prove the following.

Claim. If $\mathbf{s}_{1},\ldots,\mathbf{s}_{5}\in{\mathbb{R}}^{3}$ are points
such that

$f(\mathbf{a}_{1},\ldots,\mathbf{a}_{5},\mathbf{s}_{1},\ldots,\mathbf{s}_{5})=0\quad\text{for
all}\quad M=(\mathbf{a}_{1},\ldots,\mathbf{a}_{5})\in{\mathbb{R}}^{3\times
5},$

then $\mathbf{s}_{1}=\mathbf{s}_{2}=\cdots=\mathbf{s}_{5}$.

Having reduced the proof to the claim, we can forget about the situation of
the theorem. In principle, the claim could be proved by explicitly forming the
polynomial $f$, regarding the variables for the $\mathbf{a}_{i}$ as main
variables, extracting the coefficients (which are polynomials in the
$\mathbf{s}_{i}$-variables), and showing that the ideal generated by the
coefficients defines the variety given by
$\mathbf{s}_{1}=\mathbf{s}_{2}=\cdots=\mathbf{s}_{5}$. Unfortunately, $f$ has
15 variables and is homogeneous of degree 160. The first step towards making
the computation feasible is choosing suitable Cartesian coordinates as
follows. The vector $\mathbf{s}_{1}$ can be taken as the origin of the
coordinate system. This turns the matrix
$S:=(\mathbf{s}_{1},\ldots,\mathbf{s}_{5})\in{\mathbb{R}}^{3\times 5}$ into
$S=(\mathbf{0},\mathbf{s}_{2},\ldots,\mathbf{s}_{5})$. We can now apply QR-
decomposition, i.e., write

$S=Q\cdot\begin{pmatrix}0&b_{1}&b_{2}&b_{3}&b_{4}\\\ 0&0&b_{5}&b_{6}&b_{7}\\\
0&0&0&b_{8}&b_{9}\end{pmatrix}=:QB$ (3.2)

with $Q\in\operatorname{SO}(3)$. So using the columns of $Q$ as a new basis of
${\mathbb{R}}^{3}$, $S$ becomes the above upper triangular matrix $B$:
$(\mathbf{s}_{1},\ldots,\mathbf{s}_{5})=S=B$.

Form the matrices

$X=\begin{pmatrix}x_{1,1}&\cdots&x_{1,5}\\\ \vdots&&\vdots\\\
x_{3,1}&\cdots&x_{3,5}\end{pmatrix}\quad\text{and}\quad
Y=(Y_{i,j})=\begin{pmatrix}0&y_{1}&y_{2}&y_{3}&y_{4}\\\
0&0&y_{5}&y_{6}&y_{7}\\\ 0&0&0&y_{8}&y_{9}\end{pmatrix}$

with $x_{i,j}$ and $y_{i}$ indeterminates. With additional indeterminantes
$z_{1},\ldots,z_{5}$, form the ideal

$J:=\Bigl{(}z_{j}^{2}-\sum_{i=1}^{3}(Y_{i,j}-x_{i,j})^{2}\bigl{|}j=1,\ldots,5\Bigr{)},$

with the idea that the $z_{j}$ stand for
$\lVert\mathbf{s}_{i}-\mathbf{a}_{i}\rVert$. Modulo $J$, the product

$\prod_{\begin{subarray}{c}\varepsilon_{1},\ldots,\varepsilon_{4}=\pm 1\\\
\epsilon_{5}=1\end{subarray}}\\!\\!\\!\\!\det\left((\varepsilon_{i}z_{i}-\varepsilon_{j}z_{j})^{2}-\sum_{k=1}^{3}(x_{k,i}-x_{k,j})^{2}\right)_{i,j=1,\ldots,5}$
(3.3)

reduces to a polynomial $F(x_{1,1},\ldots,x_{3,5},y_{1},\ldots,y_{9})$ which
does not involve the $z_{j}$. If we specialize the variables in the matrix $X$
to the entries in a matrix
$M=(\mathbf{a}_{1},\ldots,\mathbf{a}_{5})=(a_{i,j})\in{\mathbb{R}}^{3\times
5}$ and the variables in $Y$ to the matrix $B$ in (3.2), we obtain

$F(a_{1,1},\ldots,a_{3,5},b_{1},\ldots,b_{9})=f(\mathbf{a}_{1},\ldots,\mathbf{a}_{5},\mathbf{s}_{1},\ldots,\mathbf{s}_{5}).$

So to prove the above claim and thus the theorem, we need to form
$F(x_{1,1},\ldots,x_{3,5},y_{1},\ldots,y_{9})$, regard the $x_{i,j}$ as main
variables and consider the ideal $L\subseteq{\mathbb{R}}[y_{1},\ldots,y_{9}]$
generated by the coefficients. Then we need to show that $L$ has
$y_{1}=\cdots=y_{9}=0$ as the only real solution, which says that all
$\mathbf{s}_{i}$ are zero and hence equal.

Alas, even after reducing the number of variables by our choice of
coordinates, computing the polynomial $F$ is still utterly impossible. What we
did instead was setting almost all of the variables $x_{i,j}$ to zero and
compute the product (3.3) with these specializations, always reducing modulo
$J$. This turns out to be possible in many cases, and extracting coefficients
gave us some generators of $L$. Doing this for many choices of specialized
variables provides ever more generators of $L$. Each time, we also reduced
modulo the generators of $L$ already known, which is permissible and
accelerates the computation. When we found a sum of squares of variables in
the ideal, we substituted this by the variables themselves since only real
solutions need to be considered.

With this technique, using random specializations of variables, we eventually
arrived at an ideal whose only solution is the origin $y_{i}=0$. We recorded
exactly which sequence of specializations led to this result and from this
produced a deterministic, reproducible procedure for verifying the claim. All
computations were done in MAGMA [?magma]. ∎

######  Remark.

Theorem 3.1 gives the theoretical justification for a procedure that detects
walls by solving the pseudo-range multilateration problem for each wall. The
following alternative method comes to mind. Since the sound traveling directly
from the loudspeaker to the microphones always arrives first, before any
echoes, the very first signals registered by the microphones must be the ones
coming directly from the loudspeaker. Therefore these can be used, without any
matching process, to determine the emission time $t$ by multilateration. Since
the echoes virtually come from the mirror points and are emitted at the same
time $t$ which is now known, the methods from ?DPWLV1 [?DPWLV1] or
?Boutin:Kemper:2019 [?Boutin:Kemper:2019] can then be used for the wall
detection.

While this alternative would presumably work in many cases, it has some
drawbacks. For one thing, obstacles may in some cases prevent one or more
microphones from hearing the direct signal from the loudspeaker. Secondly, and
perhaps more importantly, a spurious signal may be registered by one or more
of the microphones before the true signal from the loudspeaker arrives. In
this case the alternative method would mistake this spurious signal as the
direct signal. Therefore the multilateration would produce a drastically wrong
time of emission, and all subsequent wall detections, based on this erroneous
time, would become false. $\triangleleft$

## References

  * [1]
  * [2] cite.Bancroft:1985 Stephen Bancroft, _An algebraic solution of the GPS equations_ , IEEE Transactions on Aerospace and Electronic Systems 7 (1985), 56–59.
  * [3]
  * [4] cite.Beck:Pan:2012 Amir Beck, Dror Pan, _On the solution of the GPS localization and circle fitting problems_ , SIAM J. Optim. 22 (2021), 108–134.
  * [5]
  * [6] cite.magma Wieb Bosma, John J. Cannon, Catherine Playoust, _The Magma Algebra System I: The User Language_ , J. Symb. Comput. 24 (1997), 235–265.
  * [7]
  * [8] cite.Boutin:Kemper:2019 Mireille Boutin, Gregor Kemper, _A Drone Can Hear the Shape of a Room_ , SIAM J. Appl. Algebra Geometry 4 (2020), 123–140.
  * [9]
  * [10] cite.Cayley:1841 Arthur Cayley, _On a theorem in the geometry of position_ , Cambridge Mathematical Journal II (1841), 267––271.
  * [11]
  * [12] cite.Chaffee:Abel:1994 James Chaffee, Jonathan S. Abel, _On the Exact Solutions of Pseudorange Equations_ , IEEE Transactions on Aerospace and Electronic Systems 30 (1994), 1021–1030.
  * [13]
  * [14] cite.DPWLV1 Ivan Dokmanić, Reza Parhizkar, Andreas Walther, Yue M. Lu, Martin Vetterli, _Acoustic echoes reveal room shape_ , Proceedings of the National Academy of Sciences 110 (2013).
  * [15]
  * [16] cite.el2017time Youssef El Baba, Andreas Walther, Emanuël AP Habets, _Time of arrival disambiguation using the linear Radon transform_ , in: _2017 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP)_ , pp. 106–110, IEEE, 2017.
  * [17]
  * [18] cite.jager2016room Ingmar Jager, Richard Heusdens, Nikolay D Gaubitch, _Room geometry estimation from acoustic echoes using graph-based echo labeling_ , in: _2016 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP)_ , pp. 1–5, IEEE, 2016.
  * [19]
  * [20] cite.Krause:1987 Lloyd O. Krause, _A Direct Solution to GPS-Type Navigation Equations_ , IEEE Transactions on Aerospace and Electronic Systems 23 (1987), 225–232.
  * [21]
  * [22] cite.kreissig2013fast Martin Kreißig, Bin Yang, _Fast and reliable TDOA assignment in multi-source reverberant environments_ , in: _2013 IEEE International Conference on Acoustics, Speech and Signal Processing_ , pp. 355–359, IEEE, 2013.
  * [23]
  * [24] cite.li2010design Wei Li, Deng Li, Shuhui Yang, Zhiwei Xu, Wei Zhao, _Design and analysis of a new gps algorithm_ , in: _2010 IEEE 30th International Conference on Distributed Computing Systems_ , pp. 40–51, IEEE, 2010.
  * [25]
  * [26] cite.lundberg2001alternative John B Lundberg, _Alternative algorithms for the GPS static positioning solution_ , Applied mathematics and computation 119(1) (2001), 21–34.
  * [27]
  * [28] cite.sha2021reliability AZ Sha’ameri, WA Wan Aris, S Sadiah, TA Musa, _Reliability of Seismic Signal Analysis for Earthquake Epicenter Location Estimation Using 1 Hz GPS Kinematic Solution_ , Measurement 182 (2021), 109669.
  * [29]
  * [30] cite.zannini2010improved Cecilia Maria Zannini, Albenzio Cirillo, Raffaele Parisi, Aurelio Uncini, _Improved TDOA disambiguation techniques for sound source localization in reverberant environments_ , in: _Proceedings of 2010 IEEE International Symposium on Circuits and Systems_ , pp. 2666–2669, IEEE, 2010\.
  * [31]