Given two points on a unit sphere, how to express their angular difference in spherical coordinates? Given two points on a unit sphere in spherical coordinates representation ($\theta$ is longitude, $\phi$ is latitude in $[0,\pi]$), I want to express their angular difference in spherical coordinates.
In detail: given two points $p_1,p_2$ on the meridian, the first is under the second. Making the first the new north pole but leaving the meridian at $\theta =0°$ should give me a spherical coordinate with $\phi = |\phi_1 - \phi_2| $ which seems obvious, but the new longitude $\theta$ needs to become $180°$ although $\theta_1 = \theta_2 = 0$ (meridian).
If the answer is to convert both coordinates to Euclidian (X,Y,Z), then subtract and then reconvert to spherical coordinates I would not be surprised. But I would be very glad, if you could state why such a difference as in $ \mathbb{R}^{3}$ is not possible in spherical coordinates. What property do they miss?
Edit: these images should clarify my intentions.


The left image shows an angular difference. But I understand that a coordinate-wise difference is not sufficient. To show an example, please look at the right image.
$A$ becomes the new north pole. Thereby $B's$ polar coordinate changes to $\theta(B-A)$. But additionally the azimuthal angle of $B'$ becomes $\phi(B)=$ 180°. Imagine another extreme case where $\phi(B)$ started with 90°. It would end up at the same 90° after making $A$ the new north pole.   
 A: Given two points on a sphere $(x_1,y_1,z_1)=(\sin \phi_1 \cos \theta_1,\sin \phi_1 \sin \theta_1,\cos \phi_1 )$ and $(x_2,y_2,z_2)=(\sin \phi_2 \cos \theta_2,\sin \phi_2 \sin \theta_2,\cos \phi_2 )$. Take the dot product of these to obtain the cosine of the angle $\alpha$ between them relative to the center of the sphere. We have 
\begin{eqnarray*}
\cos \alpha = \cos \phi_1 \cos \phi_2 +\sin \phi_1 \sin \phi_2 ( \cos \theta_1 \cos \theta_2 +\sin \theta_1  \sin \theta_2) \\
\cos \alpha = \cos \phi_1 \cos \phi_2 +\sin \phi_1 \sin \phi_2 \cos (\theta_1 - \theta_2 ).
\end{eqnarray*}
\begin{eqnarray*}
 \alpha = \color{blue}{\cos^{-1} \left( \cos \phi_1 \cos \phi_2 +\sin \phi_1 \sin \phi_2 \cos (\theta_1 - \theta_2 )\right)}.
\end{eqnarray*}
A: The case you described is a change from one spherical coordinates system to another by exchanging the former zenith axis with a new direction expressed as a point on the sphere.
You might want to have a look at the solution of these questions: 


*

*Converting between spherical coordinate systems

*Transforming from one spherical coordinate system to another
Another method that worked for me is to do two rotations:



*

*Rotate point $B$ on the sphere around the zenith axis for $\phi(\,A\,)$ degrees. This positions $A$ on the meridian and makes $\phi(\,B\,)=\phi(\,B-A\,)$. This is a normal difference operation.

*Rotate $B$ on the sphere along the meridian around the $(\theta=90°,\phi=90°)$ axis for $\theta(\,A\,)$ degrees. This rotates $A$ along the meridian and puts it in place of the north pole. The relation between $A$ and $B$ will be preserved and your example cases are all covered. I cannot offer a formula for that, maybe others can.


If you need to implement the rotation of a point around an arbitrary axis, it is necessary to use quaternions which should be provided by your language's libraries. For C/C++ I recommend glm. I computed the rotation in XYZ-coordinates and transformed the resulting unit vector back to spherical coordinates.
A: The original question says "θ is longitude, ϕ is latitude in [0,π]", but in the associated diagram, θ is shown as colatitude, and ϕ as longitude.
Which is correct for the following discussion?
