Project 2

GPS, Conditioning, and Nonlinear Least Squares

Group Members: Patrick Bishop, Paul Mcnulty, Luke Bouck

The purpose of this project was to calculate the (x,y,z) coordinates of an object on Earth using coordinates from 4 different satellites in orbit. Since the equation for distance is not accurate for satellites, a fourth term is added to the distance equation to account for the time delay. This causes the equation to become nonlinear so solving these problems require certain methods such as multivariable newton's method. The distance equations are: \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*} where \( (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.

Problem 1

For this problem, we were given four sets of coordinates provided by four different satellites (with position measured in kilometers and time measured in seconds) and used multi-variable newton's method to solve the nonlinear system of equations given by (1). The code we used can be found here: Matlab Code .

We were given 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). We next computed the Jacobian for the system of equations given in (1) and then used this to solve multi-variable newton's method to yield the coordinates of the object on Earth to be $$ (-41.77271 km, -16.78919 km, 6370.0596 km) $$ and time delay of $$-3.201566 * 10^{-3} seconds.$$

Problem 2

For part 2, we solved for the same situation as in problem 1 but with a different method. The method we used for this problem started by subtracting each of the bottom three equations from (1) from the top equation from (1). This gave us a vector equation, $$\vec{u_x}x+\vec{u_y}y+\vec{u_z}z+\vec{w}=0.$$ We can then create an equation of x in terms of d by doing $0=\det[\vec{u_y}|\vec{u_z}|\vec{u_x}x+\vec{u_y}y+\vec{u_z}z+\vec{w}]$. The same process can be done for y in terms of d and z in terms of d. We can then substitute these equations for x, y, and z in the equation $$\sqrt{(x-A_1)^2+(y-B_1)^2+(z-C_1)^2}-c(t_1-d)=0$$, which is just a modification of the first equation from (1). We can then use bisection method to numerically solve for d and solve for x, y, and z.

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})$.

Problem 4

For problem 4, we converted the GPS coordinates given by the satellites into polar coordinates with the formulas, \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