Finding the Fall Point with Imperfect Bearings

Radio direction finding gives you a bearing—an angle from your station toward the signal source—but a single bearing only tells you the meteor fell somewhere along that line. Three stations give you three lines that should intersect at the fall point. In practice, they never quite meet. Measurement error, ionospheric refraction, and timing uncertainty mean your three lines form a small triangle of confusion rather than a point.

The surveying solution is least squares bearing intersection: find the point that minimizes the sum of squared perpendicular distances to all bearing lines. The math collapses to solving a 2×2 linear system—Ax = b, where A accumulates the weighted normal vectors from each bearing line.

In R, the matrix operations are terse:

# Stations at known positions, each reports bearing (degrees from north)
stations <- matrix(c(0,0, 10,0, 5,8.66), ncol=2, byrow=TRUE)
bearings <- c(45, 135, 270)  # degrees

A <- matrix(0, 2, 2); b <- c(0, 0)
for(i in 1:nrow(stations)) {
  theta <- bearings[i] * pi/180
  n <- c(-sin(theta), cos(theta))  # normal to bearing line
  A <- A + outer(n, n)
  b <- b + c(crossprod(n, stations[i,])) * n
}
solve(A, b)  # [5, 5] — the fall point estimate

Pascal’s 1980s surveying software would unroll the loop and use explicit arithmetic, readable but verbose:

program BearingFix;
var A11, A12, A22, b1, b2, det, x, y: real;

procedure AddBearing(sx, sy, bearing: real);
var theta, nx, ny, dot: real;
begin
  theta := bearing * pi / 180;
  nx := -sin(theta); ny := cos(theta);
  dot := nx*sx + ny*sy;
  A11 := A11 + nx*nx; A12 := A12 + nx*ny; A22 := A22 + ny*ny;
  b1 := b1 + dot*nx; b2 := b2 + dot*ny;
end;

begin
  A11 := 0; A12 := 0; A22 := 0; b1 := 0; b2 := 0;
  AddBearing(0, 0, 45);
  AddBearing(10, 0, 135);
  AddBearing(5, 8.66, 270);
  det := A11*A22 - A12*A12;
  x := (A22*b1 - A12*b2) / det;
  y := (A11*b2 - A12*b1) / det;
  writeln('Fall point: ', x:0:2, ', ', y:0:2);
end.

The Pascal version would typically live in a surveying package on a CP/M machine, processing field notebook data. The R version runs the same math but exploits matrix primitives that didn’t exist in most 1980s languages. Both give you a position estimate even when the bearings conflict—which they always do.