The goal of this project was to create methods for using satellites to determine our location. GPS satellites send out signals and know when the time was when the signal was send and know where they were when the signal was sent out. To determine our location in (x,y,z) corrdinates, we time how long the signal from the satellite took to get to our device. We then use this information to construct a system of nonlinear equations, they're spheres in fact, and use Multivariate Newton's Method or Nonlinear Least Squares to solve for our location. Because the clock in our device would be inaccurate, we would need a fourth satellite to determine a fourth variable, which would be the time delay in the clock in our device. Thus, if we have four satellites, we will use Multivariate Newton's, and if we have more than 4 satellites, we will use Nonlinear Least Squares. For more information on what specifically goes into the project, please visit this pdf.

Part 1

The goal of the first part is to take 4 satellites and determine our location from the information from those 4 satellites. For this, we will use multivariate newton's. The equations that we will use for multivariate newton's are below. \begin{eqnarray*} (x-A_1)^2+(y-B_1)^2+(z-C_1)^2 &=& c^2(t_1-d)^2 \\ (x-A_2)^2+(y-B_2)^2+(z-C_2)^2 &=& c^2(t_2-d)^2 \\ (x-A_3)^2+(y-B_3)^2+(z-C_3)^2 &=& c^2(t_3-d)^2 \\ (x-A_4)^2+(y-B_4)^2+(z-C_4)^2 &=& c^2(t_4-d)^2 \qquad (1) \end{eqnarray*} In the above equations, \( (A_i,B_i,C_i)\) are the spacial coordinates given by each individual satellite and \(t_i\) is the time the signal takes to reach the point on Earth. Also, $c$ is the speed of light in km/s. We will be solving for \( (x,y,z,d) \), where \( (x,y,z) \) is our location on earth and \( d\) is the time delay. To solve this, we will use Multivariant Newton's method. We must first have an initial guess, which is a 4 dimensional vector \(x_0\). Then function F, which is a vector valued function which is the 4 equations in (1) and takes the vectors \(x\) as an input. We finally have matrix J, which is the Jacobian matrix of F and is: $$J=\left[\begin{matrix}\frac{\partial f_1}{\partial x} & \frac{\partial f_1}{\partial y} & \frac{\partial f_1}{\partial z} & \frac{\partial f_1}{\partial d}\\ \frac{\partial f_2}{\partial x} & \frac{\partial f_2}{\partial y} & \frac{\partial f_2}{\partial z} & \frac{\partial f_2}{\partial d}\\ \frac{\partial f_3}{\partial x} & \frac{\partial f_3}{\partial y} & \frac{\partial f_3}{\partial z} & \frac{\partial f_3}{\partial d}\\ \frac{\partial f_4}{\partial x} & \frac{\partial f_4}{\partial y} & \frac{\partial f_4}{\partial z} & \frac{\partial f_4}{\partial d}\\ \end{matrix}\right]$$ where each entry in the matrix is a partial derivative of a component of F with respect to each variable. Then, we have the algorithm \begin{eqnarray*} J(x_n)v &=& F(x_n)\\ x_{n+1} &=& x_{n}-v \end{eqnarray*} To solve the matrix equation, we use gaussian elimination. This algorithm will then converge to the correct answer. Our code is in this file. The code takes an initial guess and 4 satellite vectors, where the entries of the satellite vectors are its (x,y,z) coordinates and the time it takes for the signal from the satellite to reach our gps reciever.

To test our code, we used an initial guess of (0, 0, 6370, 0), and were given the spatial coordinates measured by each satellite, (15600, 7540, 20140), (18760, 2750, 18610), (17610, 14630, 13480), and (19170, 610, 18390) (measured in km) and time intervals 0.07074, 0.07220, 0.07690, 0.07242 (measured in seconds). Our expected answer is \((x,y,z,d)=(-41.77271,-16.78919,6370.0596,-3.201556*10^{-3})\). From our code, we got, $(x,y,z,d)=(-41.772709570875,-16.789194106543,6370.059559223320,-3.201565830*10^{-3})$, which is what we wanted.

Part 2

We used the same test case as part 1. Our code does what is described in the preceding paragraph. The difference is that the input to the function part2 is that the input is a 4x4 matrix with the satellite vectors as row vectors. It then outputs the solution as a vector of $(x,y,z,d)$. Our function for part 2 is here and our function that is the above one variable equation is here. As before, we are looking for the output to be $(x,y,z,d)=(-41.77271,-16.78919,6370.0596,-3.201556*10^{-3})$, and our code produced $(x,y,z,d)=(-41.772709570836,-16.789194106532,6370.059559222558,-3.201565830*10^{-3})$, which is what we want.

Part 4

For part 4, a way of looking at the conditioning of the gps problem. For this problem, we will represent the satellite vectors in spherical coordinates, which are represented below. \begin{eqnarray*} A_i &=& \rho \cos(\phi_i) \cos(\theta_i)\\ B_i &=& \rho \cos(\phi_i) \sin(\theta_i)\\ C_i &=& \rho\sin(\phi_i) \end{eqnarray*} In the above equations, $(A_i,B_i,C_i)$ are the position vector components in rectangular coordinates, $\rho$ is the radius of the satellite from the center of the earth, $\phi$ is the angle of the satellite relative to the z axis, and $\theta$ is the angle of the satellite along the equator relative to the x axis. For our purposes, we are going to restrict $\rho=26570$, $0\leq \phi_i\leq \frac{\pi}{2}$, and $0\leq \theta_i\leq 2\pi$ for all $i=1,2,3,4$.

To test the conditioning of the specific problem, we need a known location, which will decide to be $(x,y,z,d)=(0,0,6370,.0001)$. We can then pick satellites that are decently spaced apart, but still within our restrictions. We then have to calculate the range of each satellite to our location and the time it takes for the signal to reach our location. The equations for that are \begin{eqnarray*} R_i &=& \sqrt{A_i^2+B_i^2+(C_i-6370)^2}\\ t_i &=& d+R_i/c \end{eqnarray*} Since the most probable error is from the satellites timing, we are basing the conditioning on changing all the satellite times by a combination of $\pm1*10^{-8}$, which $10^{-8}$ is about what the clocks on these satellites are capable of. We will then calculate the new position and compare it to the original position with an error vector $(\Delta x,\Delta y,\Delta z)$. To measure the error, we will use max position error $=||(\Delta x,\Delta y,\Delta z)||_{\infty}$, which is the infinity norm of the error vector. We will report that in meters. The other way we will measure error is with error magnification factor, which will be $$\frac{||(\Delta x,\Delta y,\Delta z)||_{\infty}}{c||(\Delta t_1, \Delta t_2, \Delta t_3, \Delta t_4)||_{\infty}}$$ The error magnification will give us a ratio of the max position error to the distance that light travels in the amount of time that we are changing it by, which is about 3 meters. We ran 8 runs with different combinations of time errors. The part 4 code is here. The part 4 code contains the positions of our satellite choices. The maximum position error and error maginfication values below.

Part 5

Part 6