Projekt: GPS og geometri – lineære og ikke-lineære ligninger#
Baseret på original køreplan til 01005 Matematik 1, DTU. Opdateret til Python og Mat1a+1b af Adam Greve Bregenhof
Indledning#
Hvordan ved din telefon præcis, hvor du er? Bag Google Maps og selvkørende biler ligger et netværk af ca. 24 satellitter og et fascinerende matematisk problem: GPS-positionering handler nemlig helt fundamentalt om at løse ligningssystemer, hvor ukendte koordinater og tidsfejl skal findes ud fra signaler der bevæger sig med lysets hastighed.
I dette projekt skal du dykke ned i den geometri og lineære algebra, der gør global navigation mulig. Vi starter med det ideelle tilfælde: Hvor skærer tre cirkler (eller kugler) hinanden? Men virkeligheden er dog mere kompleks. Satellitternes ure er ekstremt præcise atomure, mens dit ur i telefonen er billigt og upræcist. Denne lille tidsforskel skaber en fjerde ubekendt, som kræver en fjerde satellit og et mere avanceret system af ligninger.
Undervejs vil du arbejde med:
Afstandsligninger: Opstilling af ikke-lineære ligningssystemer baseret på Pythagoras i 2D og 3D.
Linearisering: Hvordan man transformerer komplekse cirkelligninger til simple lineære ligningssystemer ved hjælp af subtraktion.
Iterative metoder: Hvordan computere finder løsninger, når vi ikke kan isolere \(x\) og \(y\) direkte.
Python-implementering: Du skal bygge din egen “GPS-modtager” i kode, der kan beregne en position ud fra rå satellitdata og håndtere urfejl.
Projektet kombinerer klassisk geometri med moderne numeriske metoder.
Forberedelse#
Inden projektet anbefales det at have set på følgende stof:
Vektorer og koordinater i 2D og 3D (Mat1a, kapitel 1–2)
Lineære ligningssystemer og matrix-ligninger \(A\mathbf{x} = \mathbf{b}\) (Mat1a, kapitel 3)
Symmetrisk matrix, positiv definit matrix (Mat1b, kapitel 4)
Taylor-udvikling af funktioner af flere variable (Mat1b, kapitel 6)
Endelig kan man med fordel også læse mere om mindste kvadraters metode og normalligningerne \(A^T A \mathbf{x} = A^T \mathbf{b}\) her: https://data.math.au.dk/interactive/lt/mindstekvadrater.html. Metoden er måske allerede kendt fra lineær regression.
Projektmål#
Målet med projektet er at forstå den matematiske model bag GPS-positionering og implementere en komplet løsningsalgoritme i Python. Du vil:
Opstille ikke-lineære observationsligninger fra geometriske betragtninger.
Linearisere disse ligninger vha. Taylor-udvikling.
Implementere en iterativ løsningsmetode (Gauss-Newtons metode).
Udvide til overbestemte systemer med mindste kvadraters metode.
Note
I skal skrive en sammenhængende rapport uden nødvendigvis at besvare alle opgavenumrene nedenfor. Det forventes ikke, at alle opgaver besvares. I skal lade den endelige rapport styre af jeres interesser og ambitioner.
Opgaver markeret med (*) er valgfrie i den endelige rapport. Opgaver med (**) er særligt udfordrende og er tiltænkt interesserede.
Opsætning#
Følgende pakker skal installeres og importeres
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.patches as patches
from mpl_toolkits.mplot3d import Axes3D
from scipy.linalg import solve
from scipy.optimize import fsolve
from itertools import combinations
Nedenfor defineres de faste konstanter og de satellitdata, der bruges i resten af projektet.
# Lysets hastighed i m/s
c = 299_792_458
# Satellit-data: [X_s, Y_s, Z_s, pseudoafstand P_s] (alle i meter)
satellitter = np.array([
[ 4396623.907, -15219512.421, 21395963.449, 22745185], # Sat 4
[12508988.544, 10190642.686, 21073725.487, 20490898], # Sat 14
[15690450.142, -6024729.548, 20505506.585, 20673780], # Sat 16
[25180700.379, -7388281.420, 4051959.985, 23193809], # Sat 18
[-8441248.473, -20422306.240, 15063403.391, 26237952], # Sat 24
[-2062043.991, 17167206.349, 20373015.170, 22979470], # Sat 25
])
# Navne på satellitterne
sat_navne = [4, 14, 16, 18, 24, 25]
# Foreløbig startposition (X0, Y0, Z0)
x0 = np.array([3504300.0, 780800.0, 5252100.0])
# Kendt sand position
x_sand = np.array([3504320.6, 780753.5, 5252128.8])
print("Satellitdata indlæst.")
print(f"Startposition: X={x0[0]:.1f}, Y={x0[1]:.1f}, Z={x0[2]:.1f} m")
print(f"Sand position: X={x_sand[0]:.1f}, Y={x_sand[1]:.1f}, Z={x_sand[2]:.1f} m")
Modul 1 – Signalet: Fra transmissionstid til afstand#
Baggrund#
GPS-systemet består af ca. 24 satellitter i kredsløb om Jorden. Hver satellit udsender et radiosignal med information om sit eget tidsstempel og sin position. En GPS-modtager måler, hvor lang tid signalet er om at ankomme, og omregner denne tid til en afstand:
hvor \(c \approx 2{,}998 \times 10^8\ \text{m/s}\) er lysets hastighed og \(\Delta t\) er transmissionstiden.
Denne afstand kaldes en pseudoafstand (ikke en rigtig geometrisk afstand), fordi modtagerens ur ikke er perfekt synkroniseret med satellittens atomur. Vi vender tilbage til denne komplikation i Modul 4.
Alle beregninger foregår i koordinatsystemet WGS84 (World Geodetic System 1984): et tredimensionalt kartesisk system med origo i Jordens massemidtpunkt, \(x\)-aksen mod skæringen af Ækvatoren og Greenwichmeridianen, og \(z\)-aksen mod Nordpolen.
Opgave 1 – Tidsfejl og afstandsfejl#
Selv en lille fejl i modtagerens ur giver en stor fejl i positionsbestemmelsen.
a) En urfejl på \(\Delta t = 1\ \text{ms}\) giver en afstandsfejl på ca. 300 km. Bekræft dette ved at beregne \(c \cdot \Delta t\) for \(\Delta t = 10^{-3}\ \text{s}\).
b) Beregn den maksimale tilladte urfejl \(\Delta t_{\max}\), hvis vi ønsker en positionsfejl på under 10 meter.
c) Et almindeligt kvartsur har en nøjagtighed på ca. \(\pm 15\ \text{ms/dag}\). Beregn, hvor stor en positionsfejl dette svarer til per dag. Forklar ud fra din beregning, hvorfor urfejlen altid behandles som en ubekendt i GPS-systemet.
d) WGS84-koordinatsystemet er et højredrejet kartesisk koordinatsystem. Vi ved, at:
\(x\)-aksen peger mod skæringen af Ækvatoren og Greenwichmeridianen (0°N, 0°Ø),
\(z\)-aksen peger mod den geografiske Nordpol.
Brug definitionen på et højredrejet koordinatsystem til at bestemme retningen for \(y\)-aksen. Udregn hvilken bredde- og længdegrad \(y\)-aksen stikker ud af jordoverfladen ved, og find det tilsvarende sted på et verdenskort. Beskriv kort, hvad der geografisk befinder sig der.
Modul 2 – 2D-geometri og lineariseringstricket#
Baggrund#
Inden vi går til den fulde 3D-model, studerer vi det simplere 2D-tilfælde – uden urfejl. Her kender vi vores afstand til tre kendte punkter og ønsker at finde vores position.
Hvis en modtager befinder sig i \((x, y)\), og vi kender afstanden \(r_i\) til tre punkter \((x_i, y_i)\), gælder:
Geometrisk svarer dette til at finde skæringen af tre cirkler (hvis den findes). Ligningerne er ikke-lineære (kvadratiske i \(x\) og \(y\)), men der er et elegant trick til at linearisere dem. Sådanne algebraiske reduktioner vil vi se igen senere i projektet.
Opgave 2 – Linearisering ved subtraktion#
a) Udvid ligning 1 og ligning 2. Træk ligning 1 fra ligning 2. Vis, at de kvadratiske led \(x^2\) og \(y^2\) forsvinder, og at du får en lineær ligning i \(x\) og \(y\) på formen \(\alpha x + \beta y = \gamma\).
b) Gør det samme for ligning 1 trukket fra ligning 3. Du har nu to lineære ligninger med to ubekendte. Opstil det lineære system \(A\mathbf{x} = \mathbf{b}\) med: $\(A = \begin{pmatrix} 2(x_2-x_1) & 2(y_2-y_1) \\ 2(x_3-x_1) & 2(y_3-y_1) \end{pmatrix}, \quad \mathbf{b} = \begin{pmatrix} r_1^2-r_2^2+(x_2^2-x_1^2)+(y_2^2-y_1^2) \\ r_1^2-r_3^2+(x_3^2-x_1^2)+(y_3^2-y_1^2) \end{pmatrix}\)$
c) Implementér en Python-funktion, der løser dette system for et givet sæt af tre punkter og afstande. Test den med data fra nedenstående eksempel.
# Testdata: tre "sendere" i 2D og en modtager i (3, 4)
sendere_2d = np.array([[0.0, 0.0],
[6.0, 0.0],
[3.0, 7.0]])
modtager_sand_2d = np.array([3.0, 4.0])
# Beregn de "sande" afstande
afstande_2d = np.array([np.linalg.norm(modtager_sand_2d - s) for s in sendere_2d])
print("Afstande fra modtager til senderne:")
for i, (s, r) in enumerate(zip(sendere_2d, afstande_2d)):
print(f" r_{i+1} = {r:.4f} m (til punkt {s})")
Opgave 3 – Geometrisk fortolkning (*)#
a) Vis, at de tre linier \(L_{12}\), \(L_{13}\) og \(L_{23}\) (fundet ved parvis subtraktion af cirkelligningerne) alle skærer hinanden i ét fælles punkt \(P_C\). Bevis dette algebraisk.
b) Hvad repræsenterer linien \(L_{ij}\) geometrisk, når cirklerne \(i\) og \(j\) skærer hinanden i to punkter? Hvad sker der, når de ikke skærer hinanden?
c) ()** Vis, at \(P_C\) kan bestemmes ved at indføre hjælpevariablen \(s = x^2 + y^2\) og løse det resulterende lineære system i \((x, y, s)\).
Modul 3 – Kuglemetoden: Udvidelse til 3D#
Baggrund#
I den virkelige verden befinder GPS-satelliterne sig i tre dimensioner. En afstandsmåling til satellit \(i\) med position \((X_i, Y_i, Z_i)\) svarer geometrisk til, at modtageren \((X, Y, Z)\) befinder sig på en kugle med centrum \((X_i, Y_i, Z_i)\) og radius \(r_i\):
To kugler skærer i en cirkel, tre kugler skærer typisk i to punkter, og fire kugler giver (ideelt set) et entydigt punkt. I praksis kan man normalt forkaste det ene skæringspunkt, da det enten befinder sig inde i Jorden eller langt ude i rummet.
Lineariseringstricket fra 2D virker stadig: ved at trække ligninger fra hinanden forsvinder de kvadratiske led.
Opgave 4 – GPS i 3D uden urfejl#
a) Vis, at 3D-cirkelligningerne kan lineariseres ved parvis subtraktion, analogt til 2D-tilfældet. Opstil det lineære system \(A\mathbf{x} = \mathbf{b}\), hvor \(\mathbf{x} = (X, Y, Z)^T\).
b) Implementér en Python-funktion, der løser 3D-positioneringsproblemet uden urfejl, og test den på de sande satellitpositioner (kolonner 1-3 i satellitter-matricen), idet du erstatter pseudoafstandene med de sande geometriske afstande til den kendte modtagerposition \(x_{sand}\).
# Beregn sande geometriske afstande fra x_sand til satellitterne
geo_afstande = np.array([np.linalg.norm(x_sand - sat[:3]) for sat in satellitter])
print("Sande geometriske afstande til satelliterne (m):")
for i, (navn, r) in enumerate(zip(sat_navne, geo_afstande)):
print(f" Satellit {navn}: {r:.3f} m")
Modul 4 – Urfejlen: Den fjerde ubekendte#
Baggrund#
Hidtil har vi antaget, at afstandsmålingerne er eksakte geometriske afstande. I virkeligheden måler GPS-modtageren pseudoafstande \(P_s\): afstande beregnet ud fra transmissionstider, som er fejlbehæftede på grund af modtagerens upræcise ur.
Urfejlen \(dT\) (i sekunder) medfører en systematisk afstandsfejl \(w = c \cdot dT\) (i meter), der lægges til alle afstande:
Nu har vi fire ubekendte: \((X, Y, Z, w)\). Vi behøver mindst fire satellitter for at bestemme systemet entydigt. Det ikke-lineære ligningssystem er:
Opgave 5 – Observationsligninger og urfejlens rolle#
a) Forklar geometrisk, hvad det betyder at tilføje urfejlen \(w\). Forestil dig, at du trækker alle kuglers radier med den samme konstant \(w\) (og lader \(X, Y, Z\) variere frit). Hvornår skærer de “justerede” kugler hinanden i ét punkt?
b) Vis, at de første tre ligninger i systemet (de fire ligninger med fire ubekendte, minus urfejlen) ikke kan lineariseres med simpel subtraktion, når \(w\) er til stede. Hvad sker der med \(w\)-ledene, når du trækker ligning 2 fra ligning 1?
c) Beregn i Python forskellen mellem pseudoafstandene og de sande geometriske afstande for alle seks satellitter. Hvad er den typiske størrelse af urfejlen \(w\)? Svarer det til en bekymrende urfejl?
Modul 5 – Linearisering og Gauss-Newtons metode#
Baggrund#
Det ikke-lineære system fra Modul 4 løses ved hjælp af Gauss-Newtons metode (en variant af Newton-Raphson for mindste kvadraters problemer). Idéen er at linearisere ligningerne ved hjælp af første-ordens Taylor-udvikling i et startpunkt \(\mathbf{x}_0 = (X_0, Y_0, Z_0, w_0)^T\).
Lad \(f_i(\mathbf{x}) = \rho_i(\mathbf{x}) + w - P_i\), hvor \(\rho_i = \sqrt{(X-X_i)^2+(Y-Y_i)^2+(Z-Z_i)^2}\).
Taylor-udviklingen om \(\mathbf{x}_0\) giver: $\(f_i(\mathbf{x}_0 + \Delta\mathbf{x}) \approx f_i(\mathbf{x}_0) + \nabla f_i(\mathbf{x}_0)^T \Delta\mathbf{x}\)$
Sætter vi dette lig nul og kalder korrektionsvektoren \(\mathbf{x} = ({\bar X}, {\bar Y}, {\bar Z}, {\bar w})^T\), fås det lineære system:
Designmatricen \(A\) og højresiden \(\mathbf{b}\) er defineret som (for \(n\) satellitter):
hvor \(\rho_i^0 = \sqrt{(X_0-X_i)^2+(Y_0-Y_i)^2+(Z_0-Z_i)^2}\).
Løsningen \(\Delta\mathbf{x} = (\bar{X}, \bar{Y}, \bar{Z}, \bar{w})\) er en korrektion til startpunktet. Den iterative procedure er: $\(\mathbf{x}_{k+1} = \mathbf{x}_k + \Delta\mathbf{x}_k\)$
Opgave 6 – Designmatricen og dens fortolkning#
a) Vis ved at differentiere \(\rho_i = \sqrt{(X-X_i)^2+(Y-Y_i)^2+(Z-Z_i)^2}\), at de partielle afledte er: $\(\frac{\partial \rho_i}{\partial X} = \frac{X-X_i}{\rho_i}, \quad \frac{\partial \rho_i}{\partial Y} = \frac{Y-Y_i}{\rho_i}, \quad \frac{\partial \rho_i}{\partial Z} = \frac{Z-Z_i}{\rho_i}\)$
b) Vis, at de tre første kolonner i en række af \(A\) svarer til koordinaterne af en enhedsvektor fra satellit \(i\) mod modtageren. Hvad beskriver designmatricen \(A\) geometrisk?
c) Vis, at den fjerde kolonne i \(A\) altid består af \(-1\)’ere. Hvad sker der, hvis man vælger en anden startværdi for \(w_0\)?
Opgave 7 – Implementering og én iteration#
a) Implementér en Python-funktion byg_Ab(x0, w0, satellitter), der for et givet startpunkt \((X_0, Y_0, Z_0, w_0)\) beregner designmatricen \(A\) og højresiden \(\mathbf{b}\).
b) Løs det lineære system \(A\mathbf{x} = \mathbf{b}\) med de fire første satellitter (4, 14, 16, 18) og startpunktet fra \(x_0\) og \(w_0 = 0\). Hvad er korrektionsvektoren \(\Delta\mathbf{x}\)? Hvad er den opdaterede position?
c) Sammenlign den opdaterede position med den kendte sande position \(x_{sand}\). Hvor stor er fejlen?
def byg_Ab(x0, w0, satellitter):
"""
Bygger designmatricen A og højresiden b for Gauss-Newton-iterationen.
x0: (3,) array - nuværende positionsestimate (X0, Y0, Z0)
w0: float - nuværende estimate for urfejl (i meter)
satellitter:(n x 4) array - [Xs, Ys, Zs, Ps] for n satellitter
"""
# TILFØJ DIN KODE HER
raise NotImplementedError("Implementer funktionen byg_Ab for at fortsætte med Gauss-Newton-iterationerne.")
# Opgave 7b+c – skriv din løsning her
Opgave 8 – Iterativ løsning (Gauss-Newton)#
a) Implementér den fulde iterative Gauss-Newton procedure. Start fra \(x_0\) og \(w_0 = 0\). Gentag iterationen, indtil normen af korrektionsvektoren \(\|\Delta\mathbf{x}\|\) er under et stopkriterie (fx \(10^{-4}\) m).
b) Udskriv positionsestimatet og fejlen efter hver iteration. Hvor mange iterationer er nødvendige?
c) Undersøg konvergensens følsomhed over for startpunktet. Prøv at starte fra Jordens centrum \((0, 0, 0, 0)\). Hvad sker der?
tol = 0
max_iter = 0
def gauss_newton_gps(startpos, w_start, satellitter, tol=tol, max_iter=max_iter, verbose=True):
"""
Gauss-Newtons metode til GPS-positionering.
Returnerer: (pos, w, antal_iter, fejl_liste)
"""
# TILFØJ DIN KODE HER
raise NotImplementedError("Implementer gauss_newton_gps for at køre hele Gauss-Newton-iterationerne.")
def gauss_newton_gps(startpos, w_start, satellitter, tol=1e-4, max_iter=20, verbose=True):
"""
Gauss-Newtons metode til GPS-positionering.
Returnerer: (pos, w, antal_iter, fejl_liste)
"""
pos = startpos.copy().astype(float)
w = float(w_start)
if verbose:
print(f"{'Iter':>5} {'‖Δx‖ [m]':>12} {'Fejl [m]':>12}")
fejl_liste = []
for k in range(max_iter):
A, b = byg_Ab(pos, w, satellitter)
if A.shape[0] == 4:
delta = np.linalg.solve(A, b)
else:
delta = np.linalg.lstsq(A, b, rcond=None)[0]
pos += delta[:3]
w += delta[3]
norm_delta = np.linalg.norm(delta)
fejl = np.linalg.norm(pos - x_sand)
fejl_liste.append(fejl)
if verbose:
print(f"{k+1:>5} {norm_delta:>12.4f} {fejl:>12.4f}")
if norm_delta < tol:
break
return pos, w, k+1, fejl_liste
# Opgave 8b
print("=== Gauss-Newton med startposition x0 ===")
pos_gn, w_gn, n_iter, fejl_liste = gauss_newton_gps(x0, 0.0, sats_4)
print(f"\nEndelig position: X={pos_gn[0]:.3f}, Y={pos_gn[1]:.3f}, Z={pos_gn[2]:.3f}")
print(f"Urfejl: w = {w_gn:.3f} m → dT = {w_gn/c*1e6:.4f} µs")
# Opgave 8c – start fra jordens centrum
print("=== Gauss-Newton fra Jordens centrum (0,0,0) ===")
pos_centrum, w_centrum, n_iter_c, _ = gauss_newton_gps(
np.array([0.0, 0.0, 0.0]), 0.0, sats_4)
print(f"\nEndelig position: X={pos_centrum[0]:.3f}, Y={pos_centrum[1]:.3f}, Z={pos_centrum[2]:.3f}")
print(f"Fejl fra sand position: {np.linalg.norm(pos_centrum - x_sand):.4f} m")
Modul 6 – Overbestemte systemer: Mindste kvadraters metode#
Baggrund#
I praksis kan en GPS-modtager se 8–12 satellitter ad gangen. I stedet for kun at bruge 4 satellitter kan vi udnytte alle tilgængelige data til at reducere støjens indflydelse og øge præcisionen.
Når der er \(n > 4\) satellitter, har det lineariserede system \(A\mathbf{x} = \mathbf{b}\) flere ligninger end ubekendte og er generelt overbestemt (ingen eksakt løsning). Vi finder i stedet den løsning, der minimerer summen af kvadrerede residualer:
Det kan vises, at denne løsning opfylder normalligningerne:
Opgave 9 – Normalligningerne#
a) Lad \(\phi(\mathbf{x}) = \|\mathbf{b} - A\mathbf{x}\|^2 = (\mathbf{b} - A\mathbf{x})^T(\mathbf{b} - A\mathbf{x})\). Vis at \(\phi(\mathbf{x})\) er en kvadratisk form. Beregn gradienten \(\nabla_{\mathbf{x}} \phi\) og sæt den lig nul for at udlede normalligningerne \(A^T A \mathbf{x} = A^T \mathbf{b}\).
b) Vis, at matricen \(A^T A\) er symmetrisk. Argumentér for, at \(A^T A\) er positiv semi-definit, og at den er positiv definit (og dermed invertibel), når \(A\) har rang 4. Hvordan er matricen \(A^T A\) relateret til Hesse-matricen for den kvadratiske form \(\phi(\mathbf{x})\)? Brug denne relation til en anden ordens test at de stationære punkter for \(\phi(\mathbf{x})\) (under antagelse at \(A^T A\) er positiv definit).
c) Forklar, hvad residualvektoren \(\mathbf{R} = \mathbf{b} - A\hat{\mathbf{x}}\) fortæller om kvaliteten af positionsbestemmelsen.
Opgave 10 – GPS med alle seks satellitter#
a) Modificér Gauss-Newton-proceduren til at bruge alle seks satellitter (det overbestemte tilfælde). Brug np.linalg.lstsq til at løse det overbestemte lineære system i hvert iterationstrin.
b) Sammenlign den opnåede præcision med løsningen baseret på kun fire satellitter. Er der en forbedring?
c) Beregn residualvektoren \(\mathbf{R} = \mathbf{b} - A\hat{\mathbf{x}}\) for den endelige løsning. Hvad siger størrelsen af residualerne om måledatakvaliteten?
Opgave 11 – Spredningsanalyse (**)#
Når der er mere data end ubekendte, kan vi estimere usikkerheden på positionen via varians-kovariansmatricen:
hvor \(n\) er antallet af satellitter og \(m = 4\) er antallet af ubekendte.
a) Beregn \(\sigma_0^2\) og \(Q_x\) for løsningen med 6 satellitter.
b) Spredningerne på koordinaterne er \(\sigma_X = \sqrt{[Q_x]_{11}}\), \(\sigma_Y = \sqrt{[Q_x]_{22}}\), \(\sigma_Z = \sqrt{[Q_x]_{33}}\). Beregn disse og sammenlign med den faktiske fejl.
c) Hvad fortæller \(\sigma_0\) om kvaliteten af pseudoafstandsmålingerne?
Modul 7 – Perspektivering#
7.1 Relativitetsteorien og GPS#
En faktor, der ikke indgår i vores simple model men er afgørende i praksis, er Einsteins relativitetsteori. GPS-satellitter befinder sig ca. 20.200 km over Jordens overflade og bevæger sig med ca. 3,87 km/s. To relativistiske effekter påvirker satellitternees ure:
Speciel relativitetsteori: Satellittens hastighed får dens ur til at gå langsommere med ca. \(-7\ \mu\text{s/dag}\).
Generel relativitetsteori: Den svagere tyngdekraft i banen får uret til at gå hurtigere med ca. \(+45\ \mu\text{s/dag}\).
Den samlede effekt er \(+38\ \mu\text{s/dag}\): satellitterne ure springer frem ca. 38 mikrosekunder per dag i forhold til ure på Jordens overflade.
Opgave 12 – Relativistisk positionsfejl (*)#
a) Beregn den akkumulerede positionsfejl efter 1 time, hvis man ignorerede den relativistiske korrektion på \(+38\ \mu\text{s/dag}\).
b) Hvad ville fejlen være efter blot 1 minut? Diskutér, hvornår GPS-systemet ville have givet fuldstændig meningsløse positioner.
7.2 GDOP – Geometrisk fortyndning af præcision#
Designmatricen \(A\) beskriver ikke blot lineariseringen, men også geometrien i satellitkonstellationen. Hvis alle satellitterne er samlet i én del af himlen (dårlig geometri), er \(A\) dårligt konditioneret og løsningen upræcis.
GDOP (Geometric Dilution of Precision) er et mål for denne geometriske effekt: $\(\text{GDOP} = \sqrt{\text{tr}((A^T A)^{-1})}\)$
En lav GDOP (tæt på 1) er god; en høj GDOP (> 6) giver unøjagtige positioner.
Opgave 13 – GDOP-beregning (*)#
a) Beregn GDOP for den aktuelle satellitgeometri (4 og 6 satellitter).
b) Konstruér et kunstigt eksempel med 4 satellitter, der alle er samlet på den samme side af himlen (tæt på hinanden). Vis, at GDOP stiger.
def beregn_GDOP(A):
"""Beregner GDOP fra designmatricen A."""
# TILFØJ DIN KODE HER
raise NotImplementedError("Implementer beregn_GDOP for at evaluere geometrisk kvalitet af satellitkonfigurationen.")
Modul 8 – Ekstra undersøgelser (valgfrit)#
8.1 Alle kombinationer af fire satellitter (**)#
Opgave 14#
Brug alle \(\binom{6}{4} = 15\) kombinationer af fire satellitter ud af de seks tilgængelige. Beregn en position for hver kombination og sammenlign resultaterne.
a) Beregn gennemsnit og spredning på de 15 positionsestimater.
b) Diskutér, om denne metode er bedre eller dårligere end mindste kvadraters metode med alle seks satellitter.
c) Er der nogen kombinationer, der giver markant dårligere resultater? Undersøg sammenhængen med GDOP.
8.2 Visualisering af satellitpositioner (**)#
Opgave 15#
GPS-beregninger sker i 3D kartesiske koordinater (WGS84), men satellitpositioner visualiseres normalt i azimut-elevationsplot (himmelkort).
a) Konvertér satelliternes WGS84-positioner til lokale horisontale koordinater (azimut, elevation) set fra modtagerens position.
b) Visualisér satelliternes position på et polart plot (“sky plot”). Diskutér, hvad figuren fortæller om GDOP.
Appendiks A – GPS-satellitsystemet#
Det matematiske minimum (24 satellitter)#
Da GPS-systemet blev designet af det amerikanske forsvar, var det baseret på præcis 24 satellitter fordelt på 6 baner med 4 satellitter i hver. Dette antal er det matematiske minimum for at sikre, at en bruger overalt på Jorden altid kan “se” mindst 4 satellitter – nødvendigt for at bestemme en 3D-position og korrigere for urfejl.
Den nuværende situation#
I 2024 er der typisk 30–32 aktive GPS-satellitter. Ekstra satellitter fungerer som reserver og øger præcisionen.
Andre systemer (GNSS)#
Moderne telefoner bruger ofte alle tilgængelige systemer (GNSS):
System |
Land |
Ca. antal satellitter |
|---|---|---|
GPS |
USA |
31 |
Galileo |
EU |
28 |
GLONASS |
Rusland |
24 |
BeiDou |
Kina |
35 |
Appendiks B – WGS84 koordinatsystemet#
Alle GPS-beregninger foregår i WGS84 (World Geodetic System 1984):
Origo: Jordens massemidtpunkt
\(x\)-akse: Mod skæringen af Ækvatoren og Greenwichmeridianen
\(z\)-akse: Mod geografisk nordpol (Jordens rotationsakse)
\(y\)-akse: Vinkelret på \(x\) og \(z\) (højrehåndssystem)
Koordinaterne i WGS84 omregnes ofte til bredde- og længdegrader for præsentation, men alle beregninger sker i kartesiske \((X, Y, Z)\)-koordinater.
Litteratur og videre læsning#
GPS-bogen: K. Dueholm, M. Laurentzius, A. Jensen: GPS, Nyt Teknisk Forlag, 2005.
SIAM-artikel: G. Nord, D. Jabon, J. Nord: The global positioning system and the implicit function theorem, SIAM Rev., 40(3), 1998.
Mindste kvadrater: Se kursusmaterialet i Mat1b, kapitel 4.
Taylor-udvikling: Se kursusmaterialet i Mat1b, kapitel 6.
NumPy dokumentation: numpy.org
SciPy dokumentation: scipy.org