Angle between GPS coordinates I realize GPS Coordinates are spherical coordinates. However I know the earth is more of an ellipsoid. I need to compute with a fairly high degree of accuracy the pitch and yaw between two objects whose coordinates are expressed in Longitude,Latitude and Altitude.
By "pitch" I mean relative to true north, how many degrees does something at point A have to turn to face point B
By "yaw" I mean relative to the ground how many degrees do I have to tilt to face the object in question. Facing the ground could be 0 degrees, the horizon 90 degrees and the sky 180 degrees for example.
There are a few answers that start to approach this but nothing I can tell quite answers it. To determine pitch/yaw you have to realize its relative to something sitting on the outside of the ellipsoid at one of the two points. Simply using |A||B|cos(theta)=A dot b doesn't work because that gives the angle between them from the perspective of the center of the earth. Also that assumes you converted gps coordinates to XYZ first. (Which is done easily enough)
 A: Your general procedure would be something like


*

*Start with a global rectangular earth-centered coordinate system.

*For each of point A and B, find the "base point" on the GPS/WGS84 reference ellipsoid that is directly below/above the point. The reference ellipsoid has an equatorial radius of $R_e = 6\,378\,137.0\;\rm m$ and a polar radius of $R_p = 6\,356\,752.314\;\rm m$. What the latitude $\phi$ and longitude $\lambda$ of the geographical coordinates specify is the direction of the local zenith, that is, the surface normal of the reference ellipsoid at the base point. There's some algebra and trigonometry to be done here -- unless I've made a mistake somewhere, the coordinates of the base point work out to
$$ \bigl( R_e \cos(\lambda) \cos(\psi), R_e \sin(\lambda) \cos(\psi) , R_p \sin(\psi) \bigr) $$
where $\psi = \arctan\left( \frac{R_p}{R_e} \tan \phi \right)$. But better check for yourself that this makes sense.

*Find the actual points $A$ and $B$ in rectangular coordinates by going an appropriate distance straight upwards from the base point. Note that "upwards" means along the local normal defined by $\phi$ (not $\psi$) and $\lambda$ -- it is not straight away from the center of the earth. Subtract the resulting coordinates to find the direction of sight from A to B.

*Use a rotation matrix created from $(\phi_A, \lambda_A)$ to convert the difference vector $B-A$ into a local coordinate system at $A$ (with basis vectors pointing due east, north, up).

*Find your "pitch" and "yaw" by taking appropriate arctangents of the components of the line-of-sight in local coordinates.
Beware that the mean sea level can be several tens of meters above or below the reference ellipsoid surface in different parts of the world. Some GPS receivers can be configured to use a mathematical model of the mean sea level (a "geoid" which may either be a one-size-fits-all global model or a more precisely defined model for use in a single region/country) to display altitudes relative to that one instead of the reference ellipsoid. If that is the case for your input data, you will need to undo the correction in order to produce ellipsoid-relative vertical coordinates for used in the calculations sketched above.
