← Back to Logs

How GPS Actually Works: Triangulation Is a Lie

Try the interactive lab for this articleTake the quiz (6 questions · ~5 min)

If you've ever heard someone say "GPS works by triangulation," they were wrong. Not approximately wrong, not simplifying for convenience, but using the wrong word entirely. Triangulation measures angles. GPS measures distances. The correct term is trilateration, and the distinction matters because it changes what the system needs to measure, how many measurements it needs, and what error sources dominate the result.

GPS doesn't know the angle between you and a satellite. It knows (approximately) how far away the satellite is, derived from the travel time of a radio signal moving at the speed of light. Get enough of those distance measurements from satellites at known positions, and you can solve for your position. The catch: "enough" turns out to be four, not three, and the reason why reveals the deepest engineering challenge in the system.

Triangulation vs. Trilateration: Why the Distinction Matters

Triangulation is a technique from classical surveying. You measure angles from two known points to an unknown point, then use trigonometry to compute the unknown point's position. If you're standing in Athens and you know the precise direction to two landmarks at known positions, you can draw two lines along those angles and find your position at their intersection. This requires angle measurement instruments (theodolites, in the surveying world). GPS receivers don't have anything of the sort.

Trilateration works with distances, not angles. If you know your distance from three points in known positions, your location is the intersection of three spheres centred on those points. In two dimensions, two distance measurements give you two circles that intersect at (at most) two points; a third circle resolves the ambiguity. In three dimensions, three spheres intersect at (at most) two points, and typically one of those points is deep in space or underground, so you can discard it.

So in theory, three distance measurements should suffice. In practice, GPS needs four. The reason is the receiver's clock. Unlike the satellites, your phone does not carry a €50,000 atomic clock. Its quartz oscillator drifts by microseconds, and at the speed of light, one microsecond of timing error translates to 300 metres of range error. The fourth satellite measurement provides the extra equation needed to solve for the receiver's clock bias simultaneously with its three spatial coordinates. More on the math shortly.

The Satellite Constellation

The GPS constellation is managed by the United States Space Force (formerly the U.S. Air Force's 2nd Space Operations Squadron at Schriever Space Force Base). As of early 2026, there are 31 operational satellites in the constellation, though the baseline design calls for a minimum of 24. The extras provide redundancy and improved geometry.

The satellites orbit in six orbital planes, each inclined at 55 degrees to the equator, spaced 60 degrees apart in right ascension. Each plane contains at least four satellites. The orbital altitude is approximately 20,200 kilometres above the Earth's surface, which yields an orbital period of about 11 hours and 58 minutes, just under half a sidereal day. This is not a coincidence. It means each satellite's ground track repeats every sidereal day (23 hours 56 minutes), shifted by about 4 minutes relative to solar time. For the operators, this means the constellation geometry as seen from any point on Earth repeats daily, simplifying planning and analysis.

At 20,200 km altitude, the orbital velocity is approximately 3.87 km/s. Each satellite is visible from roughly 42% of the Earth's surface at any given time, and the constellation geometry guarantees that at least 4 satellites are visible from any point on Earth at any time, with typical visibility of 6 to 12 satellites. From Athens, you'll often see 8 to 10 GPS satellites simultaneously, plus additional satellites from other constellations.

The satellites themselves are substantial pieces of hardware. The current generation (GPS III and GPS IIIF) weigh about 2,200 kg at launch, carry solar panels generating around 4,500 watts, and have a design life of 15 years. Each one costs roughly €450 million to build and launch.

Atomic Clocks: Why Nanoseconds Matter

Each GPS satellite carries multiple atomic clocks, typically two rubidium and two cesium frequency standards, with newer GPS III satellites incorporating improved rubidium atomic frequency standards. These clocks define the fundamental precision of the system.

A cesium atomic clock works by exploiting the hyperfine transition frequency of cesium-133. When a cesium atom absorbs microwave radiation at exactly 9,192,631,770 Hz, it transitions between two specific energy states. The clock uses a feedback loop: a quartz oscillator generates a signal that's multiplied up to the microwave frequency range, this signal is passed through a beam of cesium atoms, and a detector measures how many atoms actually made the transition. The signal frequency is continuously adjusted to maximize the transition rate, locking the oscillator to the natural resonance frequency of cesium. This is, by definition, what a "second" is: 9,192,631,770 periods of that specific radiation.

Rubidium clocks use a similar principle with rubidium-87 atoms but achieve it through optical pumping rather than a thermal beam. They are smaller, lighter, and less power-hungry than cesium clocks, but slightly less stable over long periods. The GPS constellation uses both types for redundancy: rubidium for short-term stability (over seconds to hours), cesium for long-term stability (over days to months).

The clock accuracy on modern GPS satellites is on the order of 1 to 5 nanoseconds. Why does this matter? Because GPS measures distance by timing radio signals, and those signals travel at the speed of light: approximately 299,792,458 m/s. The math is straightforward:

distance_error = c × time_error
distance_error = 299,792,458 m/s × 1 × 10⁻⁹ s
distance_error ≈ 0.3 metres

One nanosecond of clock error produces roughly 30 centimetres of range error. Ten nanoseconds gives you 3 metres. One microsecond gives you 300 metres. If the satellite clocks drifted by just one millisecond (which an uncompensated quartz oscillator can do in minutes), the position error would be 300 kilometres. The entire system's accuracy depends on maintaining timing at the nanosecond level across a constellation of satellites orbiting at 20,200 km altitude.

The ground control segment maintains this accuracy through continuous monitoring. A network of monitoring stations around the world (including stations at Ascension Island, Diego Garcia, Kwajalein, Hawaii, and Colorado Springs) tracks each satellite's signal, compares the satellite's clock against the master control station's clock ensemble (a set of cesium fountain clocks and hydrogen masers that define GPS system time), and uploads corrections to each satellite twice daily. These corrections are broadcast as part of the navigation message, allowing receivers to compensate for residual clock drift.

Signal Structure: L1, L2, L5, and PRN Codes

GPS satellites broadcast on multiple frequencies in the L-band. Understanding the signal structure requires knowing three layers: the carrier frequency, the ranging code modulated onto it, and the navigation data message.

Carrier Frequencies

All GPS frequencies are derived from a fundamental frequency of 10.23 MHz, generated by the satellite's atomic clock:

  • L1: 1575.42 MHz (154 × 10.23 MHz). This is the primary civilian frequency. It carries the Coarse/Acquisition (C/A) code and the L1C signal on newer satellites.
  • L2: 1227.60 MHz (120 × 10.23 MHz). Originally military-only, now carries the L2C civil signal on Block IIR-M and later satellites.
  • L5: 1176.45 MHz (115 × 10.23 MHz). A newer signal specifically designed for safety-of-life applications (aviation, in particular). Only available on Block IIF and later satellites.

The choice of L-band frequencies is deliberate. Below about 1 GHz, cosmic background noise and terrestrial interference increase substantially. Above about 2 GHz, atmospheric attenuation (particularly from water vapour and oxygen) becomes problematic. The L-band sits in a sweet spot where atmospheric attenuation is low, ionospheric effects are manageable, and antenna sizes are reasonable.

PRN Codes and CDMA

Every GPS satellite broadcasts on the same frequencies simultaneously. The receiver distinguishes between satellites using Code Division Multiple Access (CDMA), where each satellite is assigned a unique Pseudo-Random Noise (PRN) code.

The C/A code on L1 is a Gold code with a length of 1,023 chips and a chipping rate of 1.023 Mchips/s, giving it a repetition period of exactly 1 millisecond. Gold codes are chosen because they have excellent cross-correlation properties: the correlation between any two different Gold codes in the GPS family is very low (approximately -24 dB), meaning the receiver can isolate one satellite's signal from all others despite them sharing the same frequency.

The P(Y) code is much longer: approximately 2.35 × 10¹⁴ chips at a chipping rate of 10.23 Mchips/s. Each satellite is assigned a unique one-week segment of this code, which repeats weekly. The P code is encrypted (the "Y" in P(Y)) using the W-code for military use, a process called Anti-Spoofing (AS).

How the Receiver Correlates

The receiver knows the exact PRN code sequence for each satellite. To acquire a satellite's signal, it generates a local replica of that satellite's PRN code and slides it in time against the incoming signal, computing the correlation at each offset. When the local replica aligns with the incoming code, the correlation spikes. This alignment gives the receiver two pieces of information:

  1. Which satellite: the PRN code that produced the correlation spike identifies the satellite.
  2. The code phase: the time offset at which correlation peaks, measured in fractions of a chip (about 1 microsecond per chip for C/A, corresponding to roughly 300 metres of range).

To refine the timing beyond the chip-level resolution, the receiver tracks the carrier phase. The L1 carrier wavelength is about 19 centimetres, so phase measurements can resolve position to millimetre levels in principle (this is the basis of RTK, discussed later).

The acquisition process is computationally expensive. The receiver must search a two-dimensional space: code phase (1,023 possible chip offsets for C/A) and Doppler frequency (the satellite's motion relative to the receiver shifts the carrier frequency by up to ±5 kHz). A brute-force search requires testing roughly 1,023 × 20 = 20,460 bins (assuming 500 Hz Doppler search steps). Modern receivers use hardware correlators that test many bins in parallel, or FFT-based acquisition that converts the correlation into a frequency-domain operation.

Navigation Message

Superimposed on the PRN code at 50 bits per second is the navigation data message. At 50 bps, it takes 12.5 minutes to transmit a complete message. The message contains:

  • Ephemeris data: precise orbital parameters for the transmitting satellite, valid for about 4 hours. This allows the receiver to compute the satellite's position at any given time.
  • Clock correction coefficients: polynomial coefficients (bias, drift, and drift rate) for correcting the satellite's clock relative to GPS time.
  • Almanac data: coarse orbital parameters for all satellites in the constellation. Less precise than ephemeris but sufficient for knowing roughly where each satellite is, which helps the receiver know which satellites to search for.
  • Ionospheric correction model parameters: coefficients for the Klobuchar model, allowing single-frequency receivers to estimate ionospheric delay.
  • UTC offset parameters: the difference between GPS time and UTC (which includes leap seconds, since GPS time does not observe them).

Pseudorange: The Core Measurement

The fundamental measurement a GPS receiver makes is the pseudorange. The name "pseudo" is critical: it's not the true range to the satellite because it includes the receiver's clock error.

Here's how it works. The satellite transmits its PRN code, timestamped by its atomic clock. The receiver generates a local replica of the same code using its own clock. By correlating the incoming signal with the local replica and finding the peak, the receiver determines the time offset between the satellite's transmission and the receiver's local time. Multiply by the speed of light, and you get a distance.

The problem is that the receiver's clock is not synchronised with the satellite's clock. If the receiver's clock is fast by 1 microsecond, every pseudorange measurement will be too long by 300 metres. If it's slow by 1 microsecond, every measurement will be too short by 300 metres. Critically, this error is the same for all satellites measured simultaneously, because it comes from the single receiver clock.

Formally, the pseudorange from receiver to satellite i is:

ρᵢ = |rᵢ - rᵤ| + c·δtᵤ + εᵢ

Where:

  • ρᵢ is the measured pseudorange to satellite i
  • rᵢ is the known position of satellite i (from ephemeris data)
  • rᵤ is the unknown receiver position (xᵤ, yᵤ, zᵤ)
  • c is the speed of light
  • δtᵤ is the unknown receiver clock bias
  • εᵢ captures all other errors (ionospheric delay, tropospheric delay, multipath, noise)

The unknowns are xᵤ, yᵤ, zᵤ, and δtᵤ: four unknowns. Each satellite provides one equation. You need at least four satellites. With more than four (the typical case), the system is overdetermined and the receiver uses a least-squares solution to find the best fit.

The Linearised Solution

The pseudorange equation is nonlinear because of the Euclidean distance term. In practice, receivers linearise it around an initial position estimate using a Taylor expansion. Let the initial estimate be (x₀, y₀, z₀, δt₀), and define the corrections Δx, Δy, Δz, Δt. The geometric range to satellite i from the estimate is:

R₀ᵢ = √((xᵢ - x₀)² + (yᵢ - y₀)² + (zᵢ - z₀)²)

The linearised observation equation becomes:

Δρᵢ = aᵢₓ·Δx + aᵢᵧ·Δy + aᵢᵤ·Δz + c·Δt

Where the direction cosines are:

aᵢₓ = -(xᵢ - x₀) / R₀ᵢ
aᵢᵧ = -(yᵢ - y₀) / R₀ᵢ
aᵢᵤ = -(zᵢ - z₀) / R₀ᵢ

And Δρᵢ = ρᵢ - R₀ᵢ - c·δt₀ is the residual pseudorange. In matrix form with n satellites:

Δρ = H · Δx

Where H is the n × 4 geometry matrix and Δx = [Δx, Δy, Δz, c·Δt]ᵀ. The least-squares solution is:

Δx = (Hᵀ H)⁻¹ Hᵀ Δρ

The receiver iterates: update the position estimate by Δx, recompute the geometry matrix, solve again, until the corrections converge below some threshold (typically a few centimetres). This usually takes 3 to 5 iterations.

Relativistic Corrections: Why Einstein Keeps Your Maps Working

GPS is one of the few engineering systems where both special and general relativity produce measurable, practically significant effects. If you turned off the relativistic corrections, GPS would accumulate roughly 11 kilometres of position error per day. Here's why.

Special Relativity: Time Dilation

The satellites orbit at approximately 3,870 m/s. Special relativity tells us that moving clocks run slower than stationary clocks. The time dilation factor is:

γ = 1 / √(1 - v²/c²)

For v = 3,870 m/s:

v²/c² = (3,870)² / (299,792,458)² ≈ 1.665 × 10⁻¹⁰

The fractional frequency shift is approximately -v²/(2c²) ≈ -8.33 × 10⁻¹¹. Over one day (86,400 seconds):

Δt_SR = -8.33 × 10⁻¹¹ × 86,400 ≈ -7.2 μs/day

The satellite clocks run slower by about 7.2 microseconds per day due to their orbital velocity. The negative sign means the satellite clocks tick slower than ground clocks from the perspective of the ground observer.

General Relativity: Gravitational Time Dilation

General relativity predicts that clocks in weaker gravitational fields (higher altitude) run faster. The gravitational potential at the satellite's orbital radius versus Earth's surface produces a fractional frequency offset:

Δf/f = (Φ_satellite - Φ_surface) / c²

Where Φ = -GM/r. The gravitational potential difference between the satellite altitude (orbital radius ≈ 26,571 km from Earth's centre) and Earth's surface (radius ≈ 6,371 km) gives:

Δf/f = (GM/c²) × (1/R_earth - 1/R_orbit)
Δf/f = (3.986 × 10¹⁴ / (299,792,458)²) × (1/6,371,000 - 1/26,571,000)
Δf/f ≈ +5.28 × 10⁻¹⁰

Over one day:

Δt_GR = +5.28 × 10⁻¹⁰ × 86,400 ≈ +45.6 μs/day

The satellite clocks run faster by about 45.6 microseconds per day due to the weaker gravitational field at their altitude.

Net Effect

The two effects partially cancel:

Δt_net = +45.6 - 7.2 = +38.4 μs/day

The satellite clocks gain roughly 38 microseconds per day relative to ground clocks. At the speed of light, 38 μs corresponds to:

38 × 10⁻⁶ × 299,792,458 ≈ 11,392 metres ≈ 11.4 km

Without relativistic corrections, GPS position errors would grow by about 11 kilometres per day, accumulating from the moment the corrections were disabled.

How the Correction Is Applied

The GPS system applies the relativistic correction in two ways. First, the satellite clock's fundamental frequency is adjusted before launch. Instead of oscillating at exactly 10.23 MHz, the atomic clocks are set to 10.22999999543 MHz, which is lower by the factory-applied fractional offset of 4.4647 × 10⁻¹⁰. When observed from the ground, the combined effect of the slower tick rate and the relativistic speedup produces an effective frequency of exactly 10.23 MHz as seen by a ground-based receiver.

Second, there is an additional periodic correction that accounts for the eccentricity of the satellite's orbit. GPS orbits are nearly circular but not perfectly so (eccentricity ≈ 0.01 to 0.02). As the satellite moves along its slightly elliptical orbit, its altitude and velocity both change, causing the relativistic effects to oscillate. This periodic correction is:

Δt_rel = -2(√(μ·a))/(c²) × e × sin(E)

Where μ is Earth's gravitational parameter, a is the semi-major axis, e is the orbital eccentricity, and E is the eccentric anomaly. The receiver computes this correction from the ephemeris data. The amplitude of this periodic term can be up to about 23 nanoseconds (roughly 7 metres of range) for the most eccentric GPS orbits.

Geometric Dilution of Precision (GDOP)

Even with perfect measurements, the accuracy of a GPS position solution depends on the geometric arrangement of the visible satellites. This effect is quantified by the Dilution of Precision (DOP) values.

Consider a simple two-dimensional analogy. If two satellites are nearly in the same direction from your position (say, both roughly north), the two range circles from those satellites intersect at a very shallow angle. A small error in either range measurement causes a large shift in the intersection point. If the two satellites are spread out (one north, one east), the circles intersect at nearly right angles, and the same range error produces a much smaller position shift.

DOP is formally derived from the geometry matrix H in the least-squares solution. The covariance matrix of the position solution is:

C = σ² × (Hᵀ H)⁻¹

Where σ² is the variance of the pseudorange measurements (assumed equal for all satellites). The diagonal elements of (Hᵀ H)⁻¹ give the DOP values:

(Hᵀ H)⁻¹ = | d₁₁  d₁₂  d₁₃  d₁₄ |
            | d₂₁  d₂₂  d₂₃  d₂₄ |
            | d₃₁  d₃₂  d₃₃  d₃₄ |
            | d₄₁  d₄₂  d₄₃  d₄₄ |

The various DOP metrics are:

  • GDOP (Geometric): √(d₁₁ + d₂₂ + d₃₃ + d₄₄), encompasses all position and time uncertainty
  • PDOP (Position): √(d₁₁ + d₂₂ + d₃₃), three-dimensional position uncertainty
  • HDOP (Horizontal): √(d₁₁ + d₂₂), horizontal position uncertainty only
  • VDOP (Vertical): √(d₃₃), vertical position uncertainty only
  • TDOP (Time): √(d₄₄), receiver clock uncertainty

The position error is approximately:

σ_position ≈ PDOP × σ_pseudorange

If your pseudorange measurement accuracy is 3 metres (typical for civilian C/A code) and PDOP is 2.0 (good geometry), your 3D position error is roughly 6 metres. With PDOP of 6.0 (poor geometry), the same pseudorange accuracy gives 18 metres of position error.

DOP value interpretation:

DOP Value Rating Description
1.0 Ideal Highest possible confidence
2.0-3.0 Excellent Suitable for all applications
4.0-6.0 Good Acceptable for navigation
7.0-8.0 Moderate Acceptable for less demanding applications
9.0-20.0 Fair Use with caution
> 20.0 Poor Measurements should be discarded

For good PDOP, you want satellites spread across the sky: some at high elevation, some at low elevation, distributed in azimuth. The worst case is all satellites clustered in the same part of the sky (common in urban canyons where buildings block low-elevation satellites on one side). Vertical DOP is almost always worse than horizontal DOP because satellites are only above you, never below, creating an inherently asymmetric geometry in the vertical axis. This is why GPS altitude accuracy is typically 1.5 to 3 times worse than horizontal accuracy.

Error Sources

Beyond clock and geometry effects, several error sources degrade the pseudorange measurement.

Ionospheric Delay

The ionosphere (roughly 80 to 1,000 km altitude) contains free electrons created by solar UV radiation. GPS signals, being electromagnetic waves, are slowed by these free electrons. The delay is proportional to the Total Electron Content (TEC) along the signal path and inversely proportional to the square of the carrier frequency:

Δt_iono = (40.3 × TEC) / (f² × c)

This frequency dependence is the key to mitigation. Dual-frequency receivers (measuring on both L1 and L2, or L1 and L5) can compute the ionospheric delay exactly, because the delay is different on each frequency but the geometric range is the same. Single-frequency receivers must rely on a broadcast model (the Klobuchar model), which corrects roughly 50 to 60% of the ionospheric delay.

The ionospheric delay varies dramatically. During quiet conditions at mid-latitudes, the delay adds about 2 to 5 metres of range error. During solar storms or at equatorial latitudes, it can exceed 30 metres. The ionosphere is one of the largest single error sources for civilian single-frequency receivers.

Tropospheric Delay

The troposphere (0 to about 12 km altitude) also delays GPS signals, but unlike the ionosphere, the tropospheric delay is not frequency-dependent. It has two components: the "dry" (hydrostatic) component, which accounts for about 90% of the total delay and is highly predictable from surface pressure measurements, and the "wet" component from water vapour, which is smaller but much harder to model because water vapour distribution is spatially and temporally variable.

The total tropospheric delay at zenith is roughly 2.3 metres and increases with lower satellite elevation angles (roughly as 1/sin(elevation), reaching 20+ metres at 5 degrees elevation). This is why receivers apply an elevation mask, typically 10 to 15 degrees, ignoring satellites near the horizon.

Multipath

Multipath occurs when the GPS signal reaches the receiver via a reflected path (bouncing off buildings, the ground, or other surfaces) in addition to the direct path. The reflected signal arrives slightly later, distorting the correlation peak and biasing the pseudorange measurement. Multipath errors can range from a few centimetres to tens of metres, depending on the reflective environment.

Multipath is particularly severe in urban environments. Standing in central Barcelona or near the Acropolis in Athens with tall structures nearby, your phone's GPS accuracy can degrade from 3 to 5 metres to 20 to 50 metres due to multipath from surrounding buildings. It is also the hardest error to model because it depends entirely on the local geometry of reflective surfaces.

Modern receivers mitigate multipath through antenna design (choke ring antennas that suppress low-elevation reflections), signal processing techniques (narrow correlator spacing, which narrows the correlation peak and reduces the influence of delayed reflections), and by using signals with higher chipping rates (the P(Y) code at 10.23 Mchips/s is much more resistant to multipath than the C/A code at 1.023 Mchips/s, because its shorter chip duration narrows the time window in which a reflected signal can affect the correlation).

Assisted GPS (A-GPS)

Anyone who used a standalone GPS receiver in the early 2000s remembers the cold start problem. You turn on the receiver, and it sits there for 30 to 45 seconds (sometimes minutes) before getting a fix. A-GPS dramatically reduces this time for mobile devices.

The Cold Start Problem

A GPS receiver performing a cold start knows nothing: not its approximate position, not the current time (to any useful precision), not which satellites are visible, and not their ephemeris data. It must:

  1. Search for satellite signals across the full code phase and Doppler space.
  2. Acquire at least one satellite.
  3. Demodulate the 50 bps navigation message to extract the almanac (which takes up to 12.5 minutes to receive in full).
  4. Use the almanac to narrow the search for other satellites.
  5. Collect ephemeris data for at least 4 satellites (each satellite's ephemeris takes 30 seconds to receive once the signal is acquired and demodulated).

A complete cold start, from power-on to first fix with no prior information, can take 30 seconds to over 2 minutes if signal conditions are poor.

Warm Start and Hot Start

A warm start occurs when the receiver has a recent almanac (less than a few weeks old) and a rough idea of its position and time (perhaps from a cell tower or from the last known position stored in memory). The almanac lets it predict which satellites should be visible and their approximate Doppler shifts, dramatically narrowing the acquisition search space. A warm start typically achieves a fix in 15 to 30 seconds.

A hot start occurs when the receiver has valid ephemeris data (less than about 4 hours old) and accurate time. It can immediately begin pseudorange measurements without waiting to demodulate the navigation message. A hot start can produce a fix in 1 to 5 seconds.

How A-GPS Works

A-GPS leverages the cellular network to provide the receiver with data it would otherwise need to demodulate from the satellite signals:

  • Almanac and ephemeris data: The cellular network downloads current ephemeris data from reference receivers operated by the network provider, delivered to the phone over a fast data connection rather than the 50 bps satellite link.
  • Approximate position: The phone's cell tower ID provides a rough position estimate (accurate to a few hundred metres in urban areas, a few kilometres in rural areas).
  • Accurate time: The cellular network provides time accurate to milliseconds.

With this assistance data, the receiver can skip the acquisition search and immediately begin tracking satellites and computing pseudoranges. This is why your phone gets a GPS fix in 2 to 3 seconds after you open a mapping app, while a standalone GPS unit takes much longer from a cold start.

A-GPS comes in two flavours. Mobile Station-Based (MS-Based): the phone receives assistance data and computes its own position. Mobile Station-Assisted (MS-Assisted): the phone sends raw pseudorange measurements to the network, which computes the position server-side. MS-Based is more common for consumer devices because it works even when the data connection drops after initial assistance.

Differential GPS and RTK: Centimetre Accuracy

Standard GPS gives you about 2 to 5 metres of horizontal accuracy with civilian signals. For applications like surveying, precision agriculture, construction, and autonomous vehicles, that's not good enough. Differential GPS (DGPS) and Real-Time Kinematic (RTK) positioning close the gap.

Differential GPS

The principle is simple. Place a GPS receiver at a precisely known position (a reference station). This reference station computes pseudoranges to all visible satellites and compares them against the ranges it should measure, given its known position. The differences are corrections that encode the combined effect of all spatially correlated errors (ionosphere, troposphere, satellite clock, satellite orbit errors) at that moment in time.

These corrections are broadcast to nearby rovers (the mobile receivers) via radio, internet, or satellite link. The rover applies the corrections to its own pseudorange measurements, cancelling errors that are common to both receivers. The key assumption: the reference station and the rover are close enough (within 100 to 200 km) that the errors are spatially correlated. Ionospheric conditions can vary significantly over large distances, limiting the effective range of DGPS corrections.

Standard DGPS achieves sub-metre accuracy (typically 0.3 to 1.0 metres).

RTK (Real-Time Kinematic)

RTK goes further by using carrier phase measurements instead of code-based pseudoranges. The L1 carrier wavelength is about 19 centimetres, so if you can count the number of complete carrier cycles between satellite and receiver (the "integer ambiguity") and measure the fractional phase, you get a range measurement with millimetre-level precision.

The challenge is determining the integer ambiguity: how many complete wavelengths fit between the satellite and the receiver? Getting this wrong by even one cycle introduces a 19 cm error. The process of resolving integer ambiguities is computationally intensive and uses sophisticated algorithms (LAMBDA, or Least-squares AMBiguity Decorrelation Adjustment, is the most common).

RTK uses a reference station within about 10 to 35 km that transmits its carrier phase observations and position to the rover. The rover performs double-differencing: differencing measurements between two satellites (to cancel receiver clock errors) and between the reference and rover (to cancel satellite clock errors). The remaining unknowns are the integer ambiguities and the rover position.

Once the ambiguities are resolved (which takes 30 seconds to a few minutes), RTK provides position accuracy of 1 to 2 centimetres horizontally and 2 to 3 centimetres vertically. This is the technology used by surveyors when you see them with GPS antennas on tripods. It's also what enables precision agriculture: tractors in the plains of Thessaly or Brandenburg can follow planting rows with 2 cm accuracy, pass after pass.

Network RTK

Rather than deploying your own reference station, Network RTK services use a network of permanently installed reference stations. In Europe, many countries operate such networks: HEPOS in Greece, SAPOS in Germany, the Ordnance Survey's OS Net in the UK. These networks generate area-correction parameters that allow RTK-level accuracy over a wider area, eliminating the need for a local base station.

The European CORS (Continuously Operating Reference Stations) networks form the backbone of modern surveying infrastructure. A surveyor in Athens can subscribe to HEPOS, receive corrections via mobile internet, and achieve 2 cm accuracy without setting up a base station.

Galileo vs. GPS vs. GLONASS vs. BeiDou

GPS is not the only game in town. Four global navigation satellite systems (GNSS) are fully operational as of 2026, and most modern receivers use signals from multiple systems simultaneously.

GPS (United States)

  • Satellites: 31 operational (24 baseline)
  • Orbital altitude: 20,200 km
  • Orbital planes: 6, inclined at 55°
  • Frequencies: L1 (1575.42 MHz), L2 (1227.60 MHz), L5 (1176.45 MHz)
  • Multiple access: CDMA
  • Civilian accuracy: ~2-5 m horizontal
  • Operational since: 1995 (full operational capability)

GLONASS (Russia)

  • Satellites: 24 operational
  • Orbital altitude: 19,100 km
  • Orbital planes: 3, inclined at 64.8°
  • Frequencies: L1 (~1602 MHz), L2 (~1246 MHz), L3 (1202.025 MHz)
  • Multiple access: FDMA (legacy), transitioning to CDMA on GLONASS-K
  • Civilian accuracy: ~3-7 m horizontal
  • Operational since: 1993 (full), restored 2011 after post-Soviet degradation

GLONASS historically used FDMA (Frequency Division Multiple Access) rather than CDMA, with each satellite transmitting on a slightly different frequency. This complicates receiver design and makes inter-system interoperability harder. The newer GLONASS-K satellites are adding CDMA signals at frequencies interoperable with GPS and Galileo. The higher orbital inclination (64.8° vs. GPS's 55°) gives GLONASS better coverage at high latitudes, which makes sense given Russia's geography.

Galileo (European Union)

  • Satellites: 28 operational (24 baseline + spares)
  • Orbital altitude: 23,222 km
  • Orbital planes: 3, inclined at 56°
  • Frequencies: E1 (1575.42 MHz), E5a (1176.45 MHz), E5b (1207.14 MHz), E6 (1278.75 MHz)
  • Multiple access: CDMA
  • Civilian accuracy: ~1-2 m horizontal (Open Service), <20 cm (High Accuracy Service)
  • Operational since: 2016 (initial services), full operational capability 2024

Europe built Galileo for several compelling reasons. First, sovereignty: relying on a foreign military system for critical infrastructure (aviation, timing, logistics) is a strategic vulnerability. GPS is operated by the U.S. military and can be degraded or disabled for civilian users at will (this happened during the 1999 Kosovo conflict, and "selective availability" deliberately degraded civilian accuracy until 2000). Second, Galileo was designed from scratch in the 2000s, allowing it to incorporate lessons learned from GPS and provide better accuracy with modern signal designs.

Galileo's E1 frequency is identical to GPS L1, and its E5a frequency is identical to GPS L5. This is deliberate: it allows combined GPS/Galileo receivers to use signals from both systems with shared RF front-end hardware. More satellites in the combined solution means better geometry and better accuracy.

Galileo's High Accuracy Service (HAS) is particularly notable. Operational since January 2023, it broadcasts precise orbit and clock corrections directly in the E6 signal, enabling accuracies below 20 centimetres without any external correction service. This is effectively embedded differential corrections, free to use, broadcast from the satellites themselves. No other GNSS provides this level of accuracy as a free, open service.

BeiDou (China)

  • Satellites: 45 operational (35 MEO + 10 GEO/IGSO)
  • Orbital altitude: 21,528 km (MEO), 35,786 km (GEO), 35,786 km (IGSO)
  • Orbital planes: 3 (MEO), inclined at 55°
  • Frequencies: B1I (1561.098 MHz), B1C (1575.42 MHz), B2a (1176.45 MHz), B3I (1268.52 MHz)
  • Multiple access: CDMA
  • Civilian accuracy: ~2-5 m horizontal (global), ~1 m (Asia-Pacific with augmentation)
  • Operational since: 2020 (BDS-3 global coverage)

BeiDou is unique in its hybrid architecture. In addition to the standard MEO constellation, it includes geostationary (GEO) satellites and inclined geosynchronous orbit (IGSO) satellites that provide enhanced coverage and accuracy over the Asia-Pacific region. The newer B1C and B2a frequencies are interoperable with GPS and Galileo.

Multi-GNSS

Modern receivers track all four constellations simultaneously. A phone in Berlin might be tracking 8 GPS, 7 Galileo, 5 GLONASS, and 6 BeiDou satellites at the same time, giving it 26 range measurements for a 4-unknown problem. This overdetermination dramatically improves accuracy, availability, and robustness. If one constellation has poor geometry from your location, the others fill the gaps.

GPS Spoofing and Jamming

GPS signals arrive at the Earth's surface incredibly weak. A GPS satellite transmits about 25 watts from 20,200 km away, resulting in a received power of roughly -130 dBm (about 10⁻¹⁶ watts). For comparison, a mobile phone transmitter puts out about 1 watt, which is 10¹³ times stronger than the received GPS signal. This extreme weakness makes GPS vulnerable to both jamming and spoofing.

Jamming

GPS jamming is conceptually simple: broadcast noise or a strong signal on the GPS frequencies to drown out the legitimate satellite signals. A device the size of a cigarette lighter, drawing a few watts, can deny GPS reception over a radius of several hundred metres. More powerful jammers can affect areas spanning kilometres.

Jamming is illegal in most jurisdictions but disturbingly common. Truck drivers use personal GPS jammers to defeat fleet tracking systems. In January 2022, a truck driver at Newark Liberty International Airport in the US was traced as the source of recurring GPS interference that affected the airport's ground-based augmentation system. The jammer was in the cab of a truck on a nearby highway.

State-level GPS jamming is routine in conflict zones. Russia has deployed GPS jamming extensively, affecting civilian aviation across the Eastern Mediterranean and the Baltic region. Pilots flying into Tallinn, Helsinki, and other northern European airports have reported GPS failures attributed to Russian jamming operations. In 2024, airlines operating in the Eastern Mediterranean, including routes into Athens, reported significant GPS interference.

Spoofing

Spoofing is more sophisticated and more dangerous than jamming. Instead of denying GPS, a spoofer broadcasts counterfeit GPS signals that mimic legitimate satellites but contain false timing or position data. A receiver that locks onto the spoofed signals will compute an incorrect position or time without any indication that anything is wrong.

A basic spoofer works by generating GPS-like signals at a higher power level than the real satellites. The spoofer gradually increases its power to capture the receiver's tracking loops, then slowly shifts the timing to drag the receiver's computed position away from its true location.

Real-world incidents:

  • 2017, Black Sea: Ships in the Black Sea near the Russian port of Novorossiysk reported GPS positions placing them at Gelendzhik airport, roughly 30 km inland. Over 20 vessels were affected simultaneously, consistent with a regional GPS spoofing operation.
  • 2019-2020, Iran: Iran spoofed the GPS of a US RQ-170 drone (the 2011 incident), and has repeatedly spoofed GPS in the Persian Gulf, causing commercial ships to report positions kilometres from their actual locations.
  • 2023-2024, Eastern Mediterranean: Widespread GPS spoofing affecting commercial aviation. Aircraft approaching airports in Beirut, Cairo, and across the region reported sudden GPS jumps of hundreds of kilometres. Some aircraft were spoofed to positions in Cairo while physically being near Beirut.

Detection Methods

Detecting spoofing is an active area of research. Methods include:

  • Multi-constellation cross-checking: A spoofer targeting GPS would need to simultaneously spoof Galileo, GLONASS, and BeiDou to fool a multi-GNSS receiver. Inconsistencies between constellations indicate spoofing.
  • Receiver Autonomous Integrity Monitoring (RAIM): By using more satellites than the minimum needed and checking for consistency, receivers can detect and exclude anomalous measurements.
  • Inertial navigation cross-checking: Comparing GPS-derived position changes with an inertial measurement unit (IMU). A sudden GPS position jump that isn't reflected in the accelerometer data indicates spoofing or signal anomaly.
  • Carrier phase monitoring: Sudden jumps in carrier phase that are inconsistent with receiver dynamics.
  • Signal power monitoring: Spoofed signals typically arrive at higher power than expected for authentic satellite signals. A receiver that monitors absolute signal power can flag anomalies.
  • Antenna-based techniques: Multi-antenna receivers can determine the direction of arrival of GPS signals. Authentic signals come from the sky; spoofed signals typically come from a terrestrial source at a different angle.

Galileo's Open Service Navigation Message Authentication (OSNMA), which became operational in 2023, adds cryptographic authentication to the navigation message. This doesn't prevent all spoofing (you can still replay signals or spoof the ranging codes), but it makes it significantly harder to inject false ephemeris or time data.

Code: Pseudorange Calculation and Position Estimation

To make the theory concrete, here are Python implementations. First, the pseudorange calculation and position estimation using iterative least squares.

import numpy as np
 
# WGS84 constants
C = 299_792_458.0  # speed of light, m/s
MU = 3.986004418e14  # Earth gravitational parameter, m³/s²
OMEGA_E = 7.2921151467e-5  # Earth rotation rate, rad/s
A_WGS84 = 6_378_137.0  # semi-major axis, m
F_WGS84 = 1 / 298.257223563  # flattening
 
def ecef_to_lla(x, y, z):
    """Convert ECEF coordinates to latitude, longitude, altitude (WGS84)."""
    b = A_WGS84 * (1 - F_WGS84)
    e2 = 1 - (b / A_WGS84) ** 2
    ep2 = (A_WGS84 / b) ** 2 - 1
    p = np.sqrt(x**2 + y**2)
    theta = np.arctan2(z * A_WGS84, p * b)
    lon = np.arctan2(y, x)
    lat = np.arctan2(
        z + ep2 * b * np.sin(theta)**3,
        p - e2 * A_WGS84 * np.cos(theta)**3
    )
    N = A_WGS84 / np.sqrt(1 - e2 * np.sin(lat)**2)
    alt = p / np.cos(lat) - N
    return np.degrees(lat), np.degrees(lon), alt
 
def compute_position(sat_positions, pseudoranges, x0=None, max_iter=10, tol=1e-4):
    """
    Estimate receiver position from satellite positions and pseudoranges.
 
    Parameters:
        sat_positions: (N, 3) array of satellite ECEF positions in metres
        pseudoranges: (N,) array of pseudorange measurements in metres
        x0: initial position estimate [x, y, z, cdt] (metres); defaults to origin
        max_iter: maximum iterations
        tol: convergence threshold in metres
 
    Returns:
        (x, y, z, cdt) estimated receiver position and clock bias (in metres)
    """
    n_sats = len(pseudoranges)
    if n_sats < 4:
        raise ValueError(f"Need at least 4 satellites, got {n_sats}")
 
    # Initial estimate: centre of Earth with zero clock bias
    if x0 is None:
        x0 = np.array([0.0, 0.0, 0.0, 0.0])
 
    x_est = x0.copy()
 
    for iteration in range(max_iter):
        # Compute predicted pseudoranges and geometry matrix
        H = np.zeros((n_sats, 4))
        delta_rho = np.zeros(n_sats)
 
        for i in range(n_sats):
            dx = sat_positions[i, 0] - x_est[0]
            dy = sat_positions[i, 1] - x_est[1]
            dz = sat_positions[i, 2] - x_est[2]
            r_pred = np.sqrt(dx**2 + dy**2 + dz**2)
 
            # Direction cosines (note the negative sign)
            H[i, 0] = -dx / r_pred
            H[i, 1] = -dy / r_pred
            H[i, 2] = -dz / r_pred
            H[i, 3] = 1.0  # clock bias coefficient
 
            # Residual: measured minus predicted pseudorange
            delta_rho[i] = pseudoranges[i] - (r_pred + x_est[3])
 
        # Least squares solution
        Ht_H = H.T @ H
        Ht_dr = H.T @ delta_rho
        delta_x = np.linalg.solve(Ht_H, Ht_dr)
 
        x_est += delta_x
 
        # Check convergence
        pos_correction = np.linalg.norm(delta_x[:3])
        if pos_correction < tol:
            break
 
    return x_est
 
def compute_dop(sat_positions, receiver_pos):
    """Compute DOP values from satellite geometry."""
    n_sats = len(sat_positions)
    H = np.zeros((n_sats, 4))
 
    for i in range(n_sats):
        dx = sat_positions[i, 0] - receiver_pos[0]
        dy = sat_positions[i, 1] - receiver_pos[1]
        dz = sat_positions[i, 2] - receiver_pos[2]
        r = np.sqrt(dx**2 + dy**2 + dz**2)
        H[i, 0] = -dx / r
        H[i, 1] = -dy / r
        H[i, 2] = -dz / r
        H[i, 3] = 1.0
 
    Q = np.linalg.inv(H.T @ H)
 
    gdop = np.sqrt(np.trace(Q))
    pdop = np.sqrt(Q[0, 0] + Q[1, 1] + Q[2, 2])
    hdop = np.sqrt(Q[0, 0] + Q[1, 1])
    vdop = np.sqrt(Q[2, 2])
    tdop = np.sqrt(Q[3, 3])
 
    return {"GDOP": gdop, "PDOP": pdop, "HDOP": hdop, "VDOP": vdop, "TDOP": tdop}
 
 
# Example: Simulated fix near Athens
 
# Receiver true position: Athens, Greece (approx ECEF)
# Lat: 37.9838°N, Lon: 23.7275°E, Alt: 100m
true_lat, true_lon, true_alt = 37.9838, 23.7275, 100.0
lat_r, lon_r = np.radians(true_lat), np.radians(true_lon)
N = A_WGS84 / np.sqrt(1 - (1 - (1 - F_WGS84)**2) * np.sin(lat_r)**2)
true_x = (N + true_alt) * np.cos(lat_r) * np.cos(lon_r)
true_y = (N + true_alt) * np.cos(lat_r) * np.sin(lon_r)
true_z = (N * (1 - F_WGS84)**2 + true_alt) * np.sin(lat_r)
 
# Simulated satellite positions (ECEF, metres) and true ranges
# 6 satellites at ~20,200 km altitude, spread across the sky
sat_ecef = np.array([
    [15_600_000, 7_540_000, 20_140_000],
    [18_760_000, 2_750_000, 18_610_000],
    [17_610_000, 14_630_000, 13_480_000],
    [19_170_000, 610_000, 18_390_000],
    [12_000_000, 18_000_000, 15_000_000],
    [22_000_000, 5_000_000, 12_000_000],
])
 
# Compute true ranges
true_ranges = np.sqrt(np.sum((sat_ecef - np.array([true_x, true_y, true_z]))**2, axis=1))
 
# Add receiver clock bias (e.g., receiver clock is 0.1ms fast = 29,979 m in pseudorange)
clock_bias_m = 29_979.2458  # 0.1 ms * c
 
# Add some noise (3m std dev, typical for C/A code)
np.random.seed(42)
noise = np.random.normal(0, 3.0, len(sat_ecef))
 
pseudoranges = true_ranges + clock_bias_m + noise
 
# Solve
solution = compute_position(sat_ecef, pseudoranges)
est_lat, est_lon, est_alt = ecef_to_lla(solution[0], solution[1], solution[2])
clock_bias_ns = solution[3] / C * 1e9
 
print(f"True position:      {true_lat:.6f}°N, {true_lon:.6f}°E, {true_alt:.1f} m")
print(f"Estimated position: {est_lat:.6f}°N, {est_lon:.6f}°E, {est_alt:.1f} m")
print(f"Clock bias:         {clock_bias_ns:.1f} ns ({solution[3]:.1f} m)")
 
# Compute DOP
dop = compute_dop(sat_ecef, solution[:3])
print(f"\nDilution of Precision:")
for k, v in dop.items():
    print(f"  {k}: {v:.2f}")

Running this produces a position estimate within a few metres of the true Athens position, with the clock bias recovered accurately.

Code: 3D Trilateration

Here's a standalone implementation of the pure trilateration problem (without clock bias), which demonstrates the geometric core of GPS positioning:

import numpy as np
 
def trilaterate_3d(p1, p2, p3, r1, r2, r3):
    """
    Solve the 3D trilateration problem: find the point(s) at distance
    r1, r2, r3 from points p1, p2, p3 respectively.
 
    Returns up to two solutions (the intersection of three spheres
    generally yields 0 or 2 points).
    """
    p1, p2, p3 = np.array(p1, dtype=float), np.array(p2, dtype=float), np.array(p3, dtype=float)
 
    # Transform to local coordinate system with p1 at origin
    # ex: unit vector from p1 to p2
    ex = (p2 - p1)
    d = np.linalg.norm(ex)
    ex = ex / d
 
    # ey: component of (p3 - p1) perpendicular to ex, normalised
    p3_p1 = p3 - p1
    i = np.dot(ex, p3_p1)
    ey = p3_p1 - i * ex
    j = np.linalg.norm(ey)
    ey = ey / j
 
    # ez: perpendicular to both
    ez = np.cross(ex, ey)
 
    # Solve in the local frame
    # x² + y² + z² = r1²
    # (x - d)² + y² + z² = r2²
    # (x - i)² + (y - j)² + z² = r3²
 
    x = (r1**2 - r2**2 + d**2) / (2 * d)
    y = (r1**2 - r3**2 + i**2 + j**2 - 2 * i * x) / (2 * j)
    z_sq = r1**2 - x**2 - y**2
 
    if z_sq < 0:
        if z_sq > -1e-6:  # numerical tolerance
            z_sq = 0
        else:
            return []  # no real solution
 
    z = np.sqrt(z_sq)
 
    # Transform back to original coordinate system
    base = p1 + x * ex + y * ey
    sol1 = base + z * ez
    sol2 = base - z * ez
 
    if z < 1e-10:
        return [sol1]
    return [sol1, sol2]
 
 
def select_surface_solution(solutions, earth_radius=6_371_000):
    """
    Given two trilateration solutions, return the one closer to
    Earth's surface.
    """
    if len(solutions) == 1:
        return solutions[0]
    d1 = abs(np.linalg.norm(solutions[0]) - earth_radius)
    d2 = abs(np.linalg.norm(solutions[1]) - earth_radius)
    return solutions[0] if d1 < d2 else solutions[1]
 
 
# Example: Trilateration from three satellites
 
# Three satellite positions (ECEF, metres)
sat1 = np.array([15_600_000, 7_540_000, 20_140_000])
sat2 = np.array([18_760_000, 2_750_000, 18_610_000])
sat3 = np.array([17_610_000, 14_630_000, 13_480_000])
 
# True receiver position (Athens, approximately)
receiver = np.array([4_617_090, 2_024_940, 3_895_790])
 
# True ranges
r1 = np.linalg.norm(sat1 - receiver)
r2 = np.linalg.norm(sat2 - receiver)
r3 = np.linalg.norm(sat3 - receiver)
 
print(f"True ranges: {r1/1000:.1f} km, {r2/1000:.1f} km, {r3/1000:.1f} km")
 
solutions = trilaterate_3d(sat1, sat2, sat3, r1, r2, r3)
print(f"\nFound {len(solutions)} solution(s):")
for idx, sol in enumerate(solutions):
    dist_from_centre = np.linalg.norm(sol) / 1000
    error = np.linalg.norm(sol - receiver)
    print(f"  Solution {idx + 1}: ({sol[0]:.0f}, {sol[1]:.0f}, {sol[2]:.0f})")
    print(f"    Distance from Earth centre: {dist_from_centre:.1f} km")
    print(f"    Error from true position: {error:.2f} m")
 
best = select_surface_solution(solutions)
print(f"\nSelected (near-surface) solution: ({best[0]:.0f}, {best[1]:.0f}, {best[2]:.0f})")
print(f"  Error: {np.linalg.norm(best - receiver):.2f} m")

This pure trilateration gives a perfect result with perfect ranges. In reality, the receiver clock bias means you never have perfect ranges, which is why the full GPS solution needs four satellites and the pseudorange formulation from the previous section.

The Full Error Budget

Putting all the error sources together, here's the approximate error budget for a single-frequency civilian GPS receiver in typical conditions:

Error Source Range Error (1σ)
Satellite clock (after correction) 1.0 m
Satellite orbit (after correction) 1.0 m
Ionospheric delay (Klobuchar model) 5.0 m
Tropospheric delay (model) 0.5 m
Multipath 1.0 m
Receiver noise 0.5 m
Total (RSS) ~5.3 m

For a dual-frequency receiver using Galileo HAS corrections, the budget shrinks dramatically: ionospheric delay drops to near zero (dual-frequency cancellation), orbit and clock errors drop to decimetres (HAS corrections), and the dominant remaining error is multipath, which is a local, environment-specific phenomenon.

The User Equivalent Range Error (UERE) is the RSS of all range error sources. The position error is then approximately UERE × DOP. With a 5 m UERE and a PDOP of 2.0, you get roughly 10 metres of 3D position error. In practice, with multi-constellation receivers and good sky visibility, typical performance is 2 to 5 metres horizontal.

Summary

GPS is a system where physics, engineering, and mathematics intersect at extraordinary precision. Atomic clocks maintain nanosecond timing. Relativistic effects that Einstein predicted a century ago are corrected in real time. Spread-spectrum signal processing extracts signals buried 20 dB below the noise floor. And a set of linear algebra equations, solved iteratively in your phone's processor, converts microsecond-level timing measurements into a position accurate to a few metres.

The fact that it works at all is remarkable. The fact that it works reliably, globally, continuously, and for free (at the point of use, courtesy of the U.S. taxpayer for GPS, and the European taxpayer for Galileo) is one of the great engineering achievements of the modern era. The next time your phone places a blue dot on a map of Paris or a navigation arrow on a highway outside Berlin, remember: four satellites, four equations, four unknowns, and a century of physics are behind that dot.