Measuring misalignment between two positions on sky

Question

Let's say I have the (RA,Dec) for two positions on the sky. I want to measure the position angle between these two positions. In other words, if I had the Cartesian coordinate vectors for these two positions (say [x,y,z]) then I think the position angle PA would simply be the dot product of the two vectors, so basically cos(PA). However, RA and Dec are spherical coordinates (corresponding to $\phi$ and $\theta$ respectively), not Cartesian coordinates. How do I get cos(PA) between the two positions on the sky if I only have their (RA,Dec) instead of [x,y,z]?

https://en.wikipedia.org/wiki/Position_angle

Position angle, usually abbreviated PA, is [a] convention for measuring angles on the sky in astronomy. The International Astronomical Union defines it as the angle measured relative to the north celestial pole (NCP), turning positive into the direction of the right ascension. In the standard (non-flipped) images this is a counterclockwise measure relative to the axis into the direction of positive declination.

Answer 1

Assuming you mean the angle between the meridian line through A and the great circle that goes through points A and B, then it goes something like this.

Define vectors from the origin to A and B assuming they lie on a unit sphere, such that $x_A= \cos \delta_A \cos \phi_A$, $ y_A = \cos \delta_A \sin \phi_A$ and $z_A= \sin \delta_A$, and similar for B. Here $\phi$ refers to right Ascension and $\delta$ is declination.

The great circles in question define planes that pass through the origin. A normal to the plane defined by OAB is given by the vector product $\vec{n_1}=\vec{A}\times \vec{B}$. Similarly a great circle passing through the O, A and the NCP $(0, 0, 1)$ has a normal of $\vec{n_2}= \vec{A}\times (0,0,1) = (y_A, -x_A, 0)$.

The angle you are looking for is the angle between these two normal vectors, which can be found from scalar product in the usual way. $$\cos \theta = \frac{\vec{n_1}\cdot \vec{n_2}}{|n_1||n_2|}=\frac{\cos \phi_{A}(y_{A}z_{B} - y_{B}z_A) + \sin \phi_{A}(x_Az_B - x_Bz_A)}{|\vec{A}\times \vec{B}|}$$

Answer 2

The position angle P of a body ($\alpha_1, \delta_1$) with respect to another body ($\alpha_2, \delta_2$) can be calculated from $$tan(P)={sin(\Delta\alpha)\over cos(\delta_2)tan(\delta_1)-sin(\delta_2)cos(\Delta\alpha)}$$ where $\Delta\alpha = \alpha_1-\alpha_2$. If the denominator is negative, the position angle lies in the range of 90 to 270 degrees.

Reference: Jean Meeus, Astronomical Algorithms, Second Edition,

Answer 3

First, it should be noted that position angle is not just defined by "two positions". The starting point is distinct from the ending point, resulting in a difference of $180^\circ$, depending on which one is first.

The previous answers were good, I just want to offer a different perspective/derivation. One way to define the position angle is that it's the angle from North counter-clockwise to the direction in question as measured on an orthographic projection that has your starting point as the origin (assuming the astronomy standard of "North up, East left").

The algebra behind the words begins by defining the unit vector/coordinates that define the starting point: $$\hat{n}_0 = \left[\begin{array}{c} \cos\alpha_0 \cos\delta_0 \\ \sin\alpha_0 \cos\delta_0 \\ \sin\delta_0 \end{array}\right],$$ with your ending position taking the same form with $0\rightarrow 1$. Orthographic projection on any vector is simply the process of projecting off the component of the vector in the direction of $\hat{n}_0$. The standard formula for this is $$\vec{v}_\perp = \vec{v} - \hat{n}_0 (\vec{v}\cdot\hat{n}_0).\tag1$$

In principle, you could apply Equation (1) numerically with $\vec{v}$ as the North pole, and then as $\hat{n}_1$, then the dot product between the normalized results will give you the position angle.

Doing it that way is probably a mistake, though. See, most positions in astronomy, $\hat{n}_0$ and $\hat{n}_1$, will be separated by a small angle, so Equation (1) will suffer badly from loss of significance. That is why it's a good idea to continue the derivation to derive a formula that doesn't have this problem.

Examining the structure of $\hat{n}_0$ is worth doing because it will simplify the algebra tremendously. It is what you get if you start with the $x$-direction unit vector, $\hat{x}$, rotate around the $y$-axis by $\delta_0$, and then rotate around the $z$-axis by $\alpha_0$. Viewed this way, finding the two vectors that are perpendicular to $\hat{n}_0$ that we need is pretty straightforward - just apply the same rotation matrices to $\hat{y}$ and $\hat{z}$. In other words, you can read them off from the columns of the rotation matrix \begin{align} R &= \left[\begin{array}{ccc} \cos\alpha_0 & -\sin\alpha_0 & 0 \\ \sin\alpha_0 & \cos\alpha_0 & 0 \\ 0 & 0 & 1 \end{array}\right] \left[\begin{array}{ccc} \cos\delta_0 & 0 & -\sin\delta_0 \\ 0 & 1 & 0 \\ \sin\alpha_0 & 0 & \cos\alpha_0 \end{array}\right] \\ & = \left[\begin{array}{ccc} \cos\alpha_0 \cos\delta_0 & -\sin\alpha_0 & -\cos\alpha_0 \sin\delta_0 \\ \sin\alpha_0 \cos\delta_0 & \cos\alpha_0 & -\sin\alpha_0 \sin\delta_0 \\ \sin\delta_0 & 0 & \cos\delta_0 \end{array}\right].\tag2 \end{align} I chose how to apply the signs to the sine functions in the two rotation matrices so that the first column would match $\hat{n}_0$.

Let's call the second column of (2) $\hat{E}'$, and the third column $\hat{N}'$. Notice that if we rotate a standard set of $x$-$y$ axes by 90 degrees counter-clockwise, then the $x$-axis corresponds to North and the $y$ to East. Thus, we can use the standard formula for the 2-dimensional component of a vector and its polar angle if we identify $\hat{n}_1\cdot\hat{N}'$ as the $x$-component and $\hat{n}_1\cdot\hat{E}'$ as the $y$. That formula is \begin{align} P &= \operatorname{atan2}\left(y,\,x\right)\\ & = \operatorname{atan2}\left(\sin\delta_1\cos\delta_0 - \sin\delta_0\cos\delta_1 \cos(\alpha_1-\alpha_0),\,\cos\delta_0\sin(\alpha_1-\alpha_0)\right). \tag3 \end{align}

Because $\cos\delta > 0$ for all declinations, you could make the formula look more like the standard one from textbooks. My own instinct is to avoid using $\tan\delta$ because it diverges near the poles. This formula will work for all $\alpha$ and $\delta$, as long as the two points are distinct. The only bit that remains is to fiddle with the $x$-like argument to make it behave well, numerically, when the points are near each other. To do so, use $\sin\delta_1\cos\delta_0 = \sin(\delta_1-\delta_0) + \sin\delta_0\cos\delta_1$ and $1 - \cos(\alpha_1-\alpha_0) = 2\sin^2\left(\frac{\alpha_1 - \alpha_0}{2}\right)$ to get $$P = \operatorname{atan2}\left(\sin(\delta_1-\delta_0) + 2\sin\delta_0\cos\delta_1 \sin^2\left(\frac{\alpha_1-\alpha_0}{2}\right),\,\cos\delta_0\sin(\alpha_1-\alpha_0)\right). \tag4$$

In principle, you'd want to investigate when it's best to use (3) or (4) numerically. In practice, I suspect that (4) will do better than (3), in terms of numerical accuracy, in the vast majority of cases that astronomers care about.