modelling Level 4

How GPS Works

A GPS receiver determines its position by measuring the time it takes for radio signals to arrive from satellites with known positions. This model derives the mathematical model from first principles: distance from time, trilateration from multiple satellites, and why you need four satellites instead of three. We'll solve the nonlinear system and explore how satellite geometry affects accuracy.

Prerequisites: distance formula, systems of equations, least squares, linearization

26 Feb 2026 Updated 12 Mar 2026 18 min read

1. The Question

How does a handheld GPS receiver know where it is on Earth?

Your phone reports your location to within a few meters. It does this by listening to radio signals from satellites 20,000 km overhead. Each satellite broadcasts:

Its position in space
The exact time the signal was transmitted

Your receiver measures when the signal arrives. From time and the speed of light, you can calculate distance. Multiple satellites give multiple distances, which should intersect at your location.

The mathematical question: How do we solve for position from distance measurements? And why do we need four satellites instead of three?

2. The Conceptual Model

Distance from Time

Radio signals travel at the speed of light: $c = 299,792,458$ m/s.

If a satellite transmits a signal at time $t_{\text{transmit}}$ and your receiver detects it at time $t_{\text{receive}}$, the signal traveled for:

\[\Delta t = t_{\text{receive}} - t_{\text{transmit}}\]

Distance:

\[d = c \cdot \Delta t\]

Example: If $\Delta t = 0.07$ seconds (70 milliseconds), then:

\[d = 299,792,458 \times 0.07 = 20,985,472 \text{ m} \approx 21,000 \text{ km}\]

This places you on a sphere of radius 21,000 km centered at the satellite.

Trilateration (Not Triangulation)

With three satellites, you have three spheres. Your position is where all three intersect.

Triangulation (used in surveying): Measure angles to known points, solve for position using geometry of triangles.

Trilateration (used in GPS): Measure distances to known points, solve for position as the intersection of spheres.

Difference: GPS receivers measure time (which gives distance), not angles.

3. Building the Mathematical Model

Known: Satellite Positions

GPS satellites broadcast their positions in ECEF coordinates (Earth-Centered, Earth-Fixed):

Origin at Earth’s center
Z-axis through North Pole
X-axis through (0° lat, 0° lon) — the intersection of equator and prime meridian
Y-axis completes the right-handed system (through 0° lat, 90° E lon)

Satellite $i$ has position $(x_i, y_i, z_i)$ in ECEF.

Unknown: Receiver Position

Receiver position: $(x, y, z)$ in ECEF.

Distance from receiver to satellite $i$:

\[r_i = \sqrt{(x - x_i)^2 + (y - y_i)^2 + (z - z_i)^2}\]

This is the geometric range — the straight-line distance through space.

Measurement: Pseudorange

The pseudorange $\rho_i$ is the measured distance based on signal travel time:

\[\rho_i = c \cdot \Delta t_i\]

In a perfect world:

\[r_i = \rho_i\]

Each satellite gives one equation. Three satellites → three equations → solve for three unknowns $(x, y, z)$.

But the world is not perfect.

4. The Clock Bias Problem

Why Three Satellites Aren’t Enough

Your receiver’s clock is not atomic. It’s a cheap quartz oscillator that drifts.

Let $b$ = receiver clock error (seconds). If your clock is 0.001 seconds fast, you’ll measure arrival times 0.001 seconds too early.

Effect on pseudorange:

\[\rho_i = c(\Delta t_i + b)\]

The measured travel time includes both the true travel time and your clock error.

This adds an unknown — the clock bias $b$ — to the system.

Four unknowns: $(x, y, z, b)$
Four equations required: Need signals from four satellites minimum.

Why This Matters

A 1-microsecond clock error causes a 300-meter position error:

\[\text{Error} = c \cdot b = 299,792,458 \times 0.000001 = 300 \text{ m}\]

GPS timing is synchronized to nanoseconds (billionths of a second). A 1-nanosecond error = 30 cm.

Satellite clocks are atomic (rubidium or cesium) and synchronized to GPS time. Your phone’s clock is not. So $b$ must be solved alongside position.

5. The Nonlinear System

Equations

For each satellite $i$, we have:

\[\sqrt{(x - x_i)^2 + (y - y_i)^2 + (z - z_i)^2} + cb = \rho_i\]

Where:

$(x, y, z)$ is the receiver position (unknown)
$b$ is the clock bias in seconds (unknown)
$(x_i, y_i, z_i)$ is satellite $i$’s position (known)
$\rho_i$ is the measured pseudorange (observed)
$c$ is the speed of light (known constant)

With $m$ satellites (typically 6–12 visible), we have $m$ equations and 4 unknowns. This is an overdetermined system — more equations than unknowns. We solve in a least-squares sense.

Why It’s Nonlinear

The square root makes this a nonlinear system. You can’t just write it as a matrix equation $\mathbf{Ax} = \mathbf{b}$ and invert.

Solution method: Linearize around an initial guess, then iterate.

6. Linearization and Iterative Solution

Step 1: Initial Guess

Start with a guess: $(x_0, y_0, z_0, b_0)$

A crude guess might be:

$(x_0, y_0, z_0)$ = center of Earth (0, 0, 0) or the centroid of satellite positions
$b_0 = 0$ (assume no clock bias initially)

Step 2: Compute Residuals

At the current guess, compute:

Predicted geometric range to satellite $i$:

\[r_{i,0} = \sqrt{(x_0 - x_i)^2 + (y_0 - y_i)^2 + (z_0 - z_i)^2}\]

Predicted pseudorange (including clock bias):

\[\hat{\rho}_i = r_{i,0} + c b_0\]

Residual (difference between observed and predicted):

\[f_i = \hat{\rho}_i - \rho_i = (r_{i,0} + c b_0) - \rho_i\]

If the guess is perfect, all residuals are zero. Otherwise, we need to adjust the guess.

Step 3: Linearize (Build the Jacobian)

Approximate each equation with a linear function near the guess.

The Jacobian matrix $\mathbf{H}$ has one row per satellite:

\[\mathbf{H}_i = \left[ \frac{\partial r_i}{\partial x}\bigg|_0, \; \frac{\partial r_i}{\partial y}\bigg|_0, \; \frac{\partial r_i}{\partial z}\bigg|_0, \; c \right]\]

Partial derivatives:

\[\frac{\partial r_i}{\partial x} = \frac{x_0 - x_i}{r_{i,0}}, \quad \frac{\partial r_i}{\partial y} = \frac{y_0 - y_i}{r_{i,0}}, \quad \frac{\partial r_i}{\partial z} = \frac{z_0 - z_i}{r_{i,0}}\]

Interpretation: These are the direction cosines — the components of the unit vector pointing from the satellite to the receiver.

The last column is $c$ because $\frac{\partial(cb)}{\partial b} = c$.

Step 4: Solve for Update

We want to find $\Delta \mathbf{x} = (\Delta x, \Delta y, \Delta z, \Delta b)$ such that:

\[\mathbf{H} \Delta \mathbf{x} \approx -\mathbf{f}\]

This is a linear least-squares problem. If $\mathbf{H}$ has more rows than columns (more satellites than unknowns), solve:

\[\Delta \mathbf{x} = -(\mathbf{H}^T \mathbf{H})^{-1} \mathbf{H}^T \mathbf{f}\]

Step 5: Update and Iterate

$x_1 = x_0 + \Delta x$ $y_1 = y_0 + \Delta y$ $z_1 = z_0 + \Delta z$ $b_1 = b_0 + \Delta b$

Repeat steps 2–5 until $

\Delta \mathbf{x}

$ is smaller than some threshold (e.g., 1 cm for position, 1 ns for time).

Typically converges in 3–5 iterations.

7. Worked Example by Hand (2D Simplified)

A full 3D GPS calculation is tedious by hand. Here’s a 2D analogue to show the mechanics.

Setup

Assume signal speed $c = 1$ (toy units).

Three transmitters at known positions:

$S_1 = (0, 0)$
$S_2 = (10, 0)$
$S_3 = (0, 10)$

True receiver position: $(x, y) = (2, 3)$ with clock bias $b = 0.5$.

Measured pseudoranges:

Calculate true ranges:

$r_1 = \sqrt{2^2 + 3^2} = \sqrt{13} = 3.606$
$r_2 = \sqrt{(2-10)^2 + 3^2} = \sqrt{73} = 8.544$
$r_3 = \sqrt{2^2 + (3-10)^2} = \sqrt{53} = 7.280$

Add clock bias ($c = 1$):

$\rho_1 = 3.606 + 0.5 = 4.106$
$\rho_2 = 8.544 + 0.5 = 9.044$
$\rho_3 = 7.280 + 0.5 = 7.780$

Iteration 1

Initial guess: $(x_0, y_0, b_0) = (1, 1, 0)$

Predicted ranges:

$r_{1,0} = \sqrt{1^2 + 1^2} = 1.414$
$r_{2,0} = \sqrt{(1-10)^2 + 1^2} = \sqrt{82} = 9.055$
$r_{3,0} = \sqrt{1^2 + (1-10)^2} = \sqrt{82} = 9.055$

Predicted pseudoranges (with $b_0 = 0$):

$\hat{\rho}_1 = 1.414$
$\hat{\rho}_2 = 9.055$
$\hat{\rho}_3 = 9.055$

Residuals:

$f_1 = 1.414 - 4.106 = -2.692$
$f_2 = 9.055 - 9.044 = 0.011$
$f_3 = 9.055 - 7.780 = 1.275$

Jacobian rows:

For $S_1 = (0, 0)$: $\left[ \frac{1-0}{1.414}, \frac{1-0}{1.414}, 1 \right] = [0.707, 0.707, 1]$

For $S_2 = (10, 0)$: $\left[ \frac{1-10}{9.055}, \frac{1-0}{9.055}, 1 \right] = [-0.994, 0.110, 1]$

For $S_3 = (0, 10)$: $\left[ \frac{1-0}{9.055}, \frac{1-10}{9.055}, 1 \right] = [0.110, -0.994, 1]$

\[\mathbf{H} = \begin{bmatrix} 0.707 & 0.707 & 1 \\ -0.994 & 0.110 & 1 \\ 0.110 & -0.994 & 1 \end{bmatrix}, \quad \mathbf{f} = \begin{bmatrix} -2.692 \\ 0.011 \\ 1.275 \end{bmatrix}\]

Solve $\mathbf{H} \Delta \mathbf{x} = -\mathbf{f}$ (using least squares or matrix inversion for the square case):

\[\Delta \mathbf{x} \approx \begin{bmatrix} 1.0 \\ 2.0 \\ 0.5 \end{bmatrix}\]

Updated guess: $(x_1, y_1, b_1) = (1, 1, 0) + (1.0, 2.0, 0.5) = (2.0, 3.0, 0.5)$

This is the true solution! In practice, you’d verify by computing residuals again (they should be near zero).

8. Computational Implementation

Below is an interactive GPS solver.

Number of satellites: 6 Measurement noise (meters): 5 m Clock bias (nanoseconds): 50 ns

Solution:

True position:

Estimated position:

Position error: m

True clock bias: ns

Estimated clock bias: ns

Iterations:

PDOP:

Try this:

Increase satellites: More satellites reduce error and improve PDOP
Add noise: See how measurement errors affect position accuracy
Change clock bias: Even large biases (100 ns) are solved correctly
Randomize satellites: Different geometries produce different PDOP values
Notice: The estimated position (pink) converges to true position (green) despite clock bias

Key insight: The fourth unknown (clock bias) is why GPS needs four satellites minimum. The system simultaneously solves for position and time.

9. Interpretation

Why the Clock Bias Trick Works

The receiver doesn’t need an atomic clock because:

All satellites share the same time reference (GPS time)
The clock bias affects all pseudoranges equally (adds the same error to each)
This creates a geometric constraint that can be solved

Analogy: If four people tell you distances to a landmark, but your measuring tape is miscalibrated by the same amount for all measurements, you can still figure out both:

Where the landmark is
How miscalibrated your tape is

Dilution of Precision (PDOP)

Position Dilution of Precision quantifies how satellite geometry affects accuracy.

Good geometry (low PDOP < 3):

Satellites spread out in all directions
Some high in sky, some near horizon
Measured distances intersect at steep angles

Bad geometry (high PDOP > 6):

Satellites clustered in one part of sky
All at similar elevations
Measured distances intersect at shallow angles

Formula: PDOP relates measurement error to position error:

\[\sigma_{\text{position}} = \text{PDOP} \times \sigma_{\text{range}}\]

Example: If ranging error is 5 m and PDOP = 2:

\[\sigma_{\text{position}} = 2 \times 5 = 10 \text{ m}\]

GPS Constellation Design

The 31 satellites in the GPS constellation are arranged so that:

6 orbital planes at 55° inclination
4–5 satellites per plane
At least 4 satellites visible from anywhere on Earth (usually 6–12)

This ensures good PDOP globally at all times.

10. What Could Go Wrong?

Atmospheric Delay

Radio signals slow down in the ionosphere (charged particles) and troposphere (water vapor).

Effect: Signals take longer than expected → pseudoranges are too large → position error.

Mitigation:

Dual-frequency receivers: Use two frequencies (L1 and L2) to estimate ionospheric delay
Atmospheric models: Predict tropospheric delay from temperature/pressure
Differential GPS: Use corrections from a nearby base station

Multipath

Signals bounce off buildings, terrain, or water before reaching the antenna.

Effect: Reflected signals travel farther → pseudoranges too large → position error (typically 5–10 m in urban areas).

Mitigation:

Antenna design: Reject signals from below (likely reflected)
Signal processing: Correlate to detect direct vs. reflected signals
Carrier phase tracking: More resistant to multipath than code phase

Poor Satellite Geometry

If all visible satellites are in one part of the sky (high PDOP):

Position error is amplified
Vertical accuracy suffers most (satellites overhead improve HDOP but not VDOP)

Mitigation: Wait for better satellite configuration, or use multi-GNSS (GPS + GLONASS + Galileo + BeiDou).

Relativistic Effects

Satellites are moving at 3.9 km/s in weaker gravity than Earth’s surface.

Special relativity: Moving clocks run slower (−7 μs/day)
General relativity: Clocks in weaker gravity run faster (+45 μs/day)
Net effect: Satellite clocks run fast by +38 μs/day

Without correction, this would cause 10 km/day position error!

Solution: Satellite clocks are pre-adjusted to run at 10.22999999543 MHz instead of 10.23 MHz exactly.

11. Extension: RTK and Carrier Phase

Standard GPS uses code phase (timing of pseudorandom noise pattern) → 1–5 m accuracy.

RTK (Real-Time Kinematic) uses carrier phase (actual radio wave oscillations) → 1–2 cm accuracy.

Carrier phase ambiguity:

The carrier phase measurement is:

\[\Phi_i = \frac{1}{\lambda}(r_i + cb) + N_i\]

Where:

$\lambda$ is the carrier wavelength (19 cm for GPS L1)
$N_i$ is an integer number of full wavelengths (ambiguity)

Challenge: Solve for both position and integer ambiguities simultaneously.

Method: Use a nearby base station with known position to form differential measurements that eliminate common errors, then resolve integer ambiguities using statistical methods (LAMBDA algorithm).

12. Math Refresher: Least Squares

Why Least Squares?

With more equations than unknowns ($m > 4$), the system is overdetermined. There’s no exact solution that satisfies all equations perfectly (due to measurement noise).

Least squares finds the solution that minimizes the sum of squared residuals:

\[\min_{\mathbf{x}} \sum_{i=1}^{m} f_i^2(\mathbf{x})\]

Matrix Form

For linear system $\mathbf{Hx} = \mathbf{b}$, the least-squares solution is:

\[\mathbf{x} = (\mathbf{H}^T \mathbf{H})^{-1} \mathbf{H}^T \mathbf{b}\]

This minimizes $|\mathbf{Hx} - \mathbf{b}|^2$.

Covariance: If measurement errors have variance $\sigma^2$, then:

\[\text{Cov}(\mathbf{x}) = \sigma^2 (\mathbf{H}^T \mathbf{H})^{-1}\]

The diagonal elements are the variances of $x$, $y$, $z$, and $b$.

Summary

GPS determines position by measuring signal travel time from satellites
Trilateration: Each satellite measurement places you on a sphere
Four satellites needed: 3 for position (x, y, z) + 1 for clock bias (b)
Clock bias is solved alongside position — no atomic clock required in receiver
Nonlinear system solved by linearization and iteration
PDOP quantifies how satellite geometry affects accuracy
Real GPS corrects for atmospheric delay, relativistic effects, and multipath
RTK uses carrier phase for centimeter-level accuracy

Next steps: Further exploration of this subject could include:

Differential GPS and SBAS (Satellite-Based Augmentation Systems)
Kalman filtering for continuous position tracking
Integration with INS (Inertial Navigation Systems)
Multi-GNSS constellation usage