The linear algebra behind the Global Postitioning System
Contents
Symbolic Variables
syms x y z zt s
The Input Data from Four Satalites
disp('Data received from four satelites');
SatXYZ = [
.94, 1.2, 2.3;
0, 1.35, 2.41;
2.3, .94, 1.2;
2.41, 0, 1.35;
]
SendTime = [1.23; 1.08; .74; .23]
Data received from four satelites
SatXYZ =
0.9400 1.2000 2.3000
0 1.3500 2.4100
2.3000 0.9400 1.2000
2.4100 0 1.3500
SendTime =
1.2300
1.0800
0.7400
0.2300
The non-linear expressions
for i = 1:4
NonLinearExpression{i} = (x-SatXYZ(i,1))^2 + (y-SatXYZ(i,2))^2 + ...
(z-SatXYZ(i,3))^2 - .22*(zt-SendTime(i))^2;
disp(NonLinearExpression{i})
end
(x - 47/50)^2 + (y - 6/5)^2 + (z - 23/10)^2 - (11*(zt - 123/100)^2)/50
(y - 27/20)^2 + (z - 241/100)^2 - (11*(zt - 27/25)^2)/50 + x^2
(x - 23/10)^2 + (y - 47/50)^2 + (z - 6/5)^2 - (11*(zt - 37/50)^2)/50
(x - 241/100)^2 + (z - 27/20)^2 - (11*(zt - 23/100)^2)/50 + y^2
Obtain a linear system from the non-linear one
for i = 1:3;
LinearExpression{i} = simplify(NonLinearExpression{1}-NonLinearExpression{i+1});
AugmentedMatrix(i,:) = flipud(coeffs(LinearExpression{i})');
disp(LinearExpression{i})
end
(3*y)/10 - (47*x)/25 + (11*z)/50 + (33*zt)/500 - 9323/100000
(68*x)/25 - (13*y)/25 - (11*z)/5 + (539*zt)/2500 - 106183/500000
(147*x)/50 - (12*y)/5 - (19*z)/10 + (11*zt)/25 - 1691/5000
Create the augmented matrix from the expression
AugmentedMatrix(:,end) = -AugmentedMatrix(:,end);
disp(double(AugmentedMatrix));
-1.8800 0.3000 0.2200 0.0660 0.0932
2.7200 -0.5200 -2.2000 0.2156 0.2124
2.9400 -2.4000 -1.9000 0.4400 0.3382
Find the solution to the linear system
B = double(rref(AugmentedMatrix))
xx = B(1,5) - B(1,4)*s;
yy = B(2,5) - B(2,4)*s;
zz = B(3,5) - B(3,4)*s;
B =
1.0000 0 0 -0.0782 -0.0875
0 1.0000 0 -0.1538 -0.1059
0 0 1.0000 -0.1583 -0.1797
The linear system has infinite solutions, one of which solves the non-linear system.
Find the value of s which solves the non-linear system
P = subs(NonLinearExpression{1},{x,y,z,zt},{xx,yy,zz,s});
S = double(coeffs(P));
Poly = [num2str(S(3)),'x^2 + ',num2str(S(2)),'x + ',num2str(S(1))]
s = double(solve(P))
t = max(s);
X = subs(xx,tt)
Y = subs(yy,tt)
Z = subs(zz,tt)
Poly =
-0.16519x^2 + -0.80611x + 8.5771
s =
5.1677
-10.0477
X =
0.3033
Y =
0.6630
Z =
0.6118
Convert to long/lat
r = sqrt(X^2+Y^2+Z^2)
long = 180/pi*atan(Y/X);
lat = 90 - 180/pi*acos(Z/r);
LongLat = [long,lat]
disp('In a left-handed coordinate system, the ship is located just off the coast of Rhode Island. ');
r =
0.9518
LongLat =
65.4152 40.0004
In a left-handed coordinate system, the ship is located just off the coast of Rhode Island.