GNSS as an optimization problem

Thu 20 July 2023
Satellite optimized (Photo credit: Wikipedia)

GNSS's basic functionality works as follow: Satellites sends their position and their local time frequently. The time is synchronized between satellites thanks to the most accurate time keeping devices ever made (atomic clock). Each GNSS constellation implements its own time (e.g. GPS time). The basic problem is to define the distance to each satellite. Knowing their positions, one can then deduce its position. I decided to tackle this problem by first finding the current time. This time is constraint by signal speed (let's assume it is light speed). Thus I decided to formulate this problem as an optimization problem: find the earliest time so that we can receive all messages sent by satellites. The constraints are quite easy: we must be located in light cones with its apex at the satellite when it send its message.

light cone


Satellite i signal its position \(x_i\) at \(t_i\)

Being in its light cone means being in position \(x\) at time \(t\) such that

\begin{equation*} \frac{\|x - x_i\|}{\|t - t_i\|} \leq c \end{equation*}

Moreover, causality implies \(t \geq t_i\).

Optimisation problem

The goal is to find intersection of cones, i.e. the fastest way for signals to travel (assumes mesages travel at light speed).

The solution of the following optimisation problem find both \(t\) and \(x\)

\begin{equation*} \begin{cases} \min_{x, t} & t \\ s.t. & t > t_i \qquad \forall \, t_i\\ & \|x - x_i\| \leq c (t-t_i) \qquad \forall \, t_i, x_i \\ \end{cases} \end{equation*}

Solutions may not be unique if there is not enough \((x_i, t_i)\) or if x and multiple satellites are aligned (light cone intersections are lines and not a denombrable number of points)

Python transcription

Such an optimization problem can be solved using scipy.optimize.

I decided to represent the event "the satellite send its position \(\vec{x_s}\) at time \(t_s\)" by the concatenation of \(\vec{x_s}\) and \(t_s\).

First, let's plot some light cones:

def light_cone(s, t):
    c = [[s - i, t + i] for i in range(10)]
    c.extend([[s + i, t + i] for i in range(10)])
    return (np.array(c).transpose())

s1, t1 = -5, 0
s2, t2 = 4, 1

c1 = light_cone(s1, t1)
c2 = light_cone(s2, t2)

plt.text(s1, t1, "S1")
plt.text(s2, t2, "S2")
plt.title("light cones")

Then let's dive into the optimization problem.

First, the satellites:

s1 = np.array([s1])
s2 = np.array([s2])
satelittes = [(s1, t1), (s2, t2)]

Now the cost function (quite easy):

cost = lambda s: s[-1]

The linear constraints can be written as:

\begin{equation*} \begin{pmatrix} t_0\\ \vdots\\ \\ t_n \end{pmatrix}% <% \begin{pmatrix} 0 & \cdots & 0 & 1 \\ \vdots & \ddots & \vdots & \vdots \\ 0 & \cdots & 0 & 1 \\ \end{pmatrix} \begin{pmatrix} x_0\\ \vdots\\ x_n\\ t \end{pmatrix} \end{equation*}

and in python:

from scipy.optimize import LinearConstraint

ncons = len(satelittes)
A = np.zeros((ndim + 1, ncons))
A[ndim][:] = np.ones((ncons,))
A = A.transpose()
lower = np.array([x[-1] for x in satelittes]) # [t0, ... t_n]
upper = np.array([np.inf for _ in range(ncons)])

linear_constraints = LinearConstraint(A, lb=lower, ub=upper)

and then, the non-linear constraints:

from scipy.optimize import NonLinearConstraint

c = 1
non_linears = []

def factory(s):
    fun = lambda x: (c * abs(x[-1] - s[-1]) - np.linalg.norm(s[:-1] - x[:-1]))
    return fun

for s in satelittes:
    xs = s[:-1]
    ts = s[-1]

    fun = factory(s)

Let's define an arbitrary initial value for the optimizer. The optimization is faster if this value already respects the constraints, i.e. the point is already in all light cones.

A simple way to have one is to take a point far in the future.

x0 = satelittes[0].copy()
x0[-1] = c * (100 + max(abs(s[-1]) for s in satelittes))

Now, we can optimize:

from scipy.optimize import minimize

result = minimize(cost, x0=x0, constraints=constraints)

x0, t0 = result.x[:-1], result.x[-1]

Final thought

This optimization problem is only one possibility to find a position using GNSS. I don't know how it is done in real life and there may be easier way to solve this problem.

Category: maths Tagged: python maths optimization

Differential equation in python

Sat 04 April 2020
Second order differential equation (Photo credit: Wikipedia)

In python, differential equations can be numerically solved thanks to scipy [1]. Is usage is not as intuitive as I expected.

Simple equation

Let's start small. The first equation will be really simple:

\begin{equation*} \frac{\partial{f}}{\partial{t}} = a \times f …

Category: maths Tagged: python maths equation

Read More

Zombie propagation

Sat 21 March 2020
Zombie favorite food warning (Photo credit: wikipedia)

I recently read a paper [1] trying to model a disease propagation. I wanted to play with this model.

The model

The model is know as "SIR" as it divide the population into 3 groups:

  • S: suceptible to become a zombie
  • I: infected …

Category: maths Tagged: python maths zombie

Read More

Music reflection: timeline and publications

Sat 18 January 2014

[caption id="" align="alignright" width="350"]English: artificial omni-directional sound sou... sound research material (Photo credit: Wikipedia)[/caption]

I just discover google music timeline. That leads me to several publications from google about music.

I wonder why does google work on music. Google is a commercial company whose business model is based on ads and …

Category: maths Tagged: Business model Google Play Publication search engine Timeline maths music reflections

Read More
Page 1 of 1