My friend Sam is part of a team building a Gravitational Wave Treasure Map. Due to the way gravitational wave telescopes (e.g., LIGO) work, gravitational wave events can only be roughly localized on the sky. For example, gravitational wave event S190923y could have come from anywhere in the dark red regions here:

As soon as a detection is made, alerts are sent out to the astronomical community. Teams of astronomers around the world then search the indicated region for new points of light, in hopes of finding an optical (electromagnetic) counterpart. They report back which regions of the sky they imaged, and the treasure map shows where to focus further searches.

In the course of developing this map, he ran in to a classic coordinate singularity problem, and I helped him work out the math behind it.

The singularity¶

Astronomers locate objects in angular coordinates, right ascension (RA) and declination (Dec)—analogous to longitude and latitude respectively. The problem is that, when your declination goes to $\pm 90^\circ$, your right ascension is undefined. (Just like when you're standing at the North Pole, you don't have a particular longitude.)

So, say you have a telescope that images a $10^\circ \times 10^\circ$ field of view. That's 10 degrees from the telescope pointing location, not 10 degrees RA by 10 degrees Dec. It's easy to see why: define a square centered at $(0^\circ, 0^\circ)$:

And try to move it up to the pole:

If you have centered on $(0^\circ, 85^\circ)$, you'll have a corner at $(-5^\circ, 90^\circ)$ and $(5^\circ, 90^\circ)$. But, of course, RA is undefined when Dec is $90^\circ$! In other words, the two top corners end up occupying the same point, and your footprint is clearly not $10^\circ$ by $10^\circ$ any more.

My friend encountered this when he tried to place detector footprint outlines at the pole for a visualization. There were division by zero errors, wild distortions... basically, it didn't work.

It turns out there's no simple formula to move your polygon around in angular coordinates while preserving surface area (solid angle).

Transforming into a higher dimensional space¶

Actually, that's not entirely true. There's no simple 2-D formula. If you convert into 3D Cartesian coordinates, you can apply standard rotation matrices to move your footprint where-ever you like.

If we have $\text{RA} \in [0^\circ,\ 360^\circ]$ and $\text{Dec} \in [-90^\circ,\ 90^\circ]$, we can convert to Cartesian coordinates by first defining

$$ \begin{align*} \phi &= (90^\circ - \text{Dec}) \times \frac{\pi}{180^\circ}\\ \theta &= \text{RA} \times \frac{\pi}{180^\circ}\\ \end{align*} $$

Then using the usual conventions, our coordinates are

$$ \begin{align*} x &= r \cos \theta \sin \phi \\ y &= r \sin \theta \sin \phi \\ z &= r \cos \phi \end{align*} $$

but we only care about the surface of the unit sphere, so we always have $r = 1$. Let's define the functions to convert back and forth:

With the conventions we've chosen, $(\text{RA},\ \text{Dec}) = (0^\circ,\ 0^\circ)$ should be along the X axis. Check and see:

Slight imprecisions in the conversion of the Dec coordinate led to a negligible nonzero $z$ component, but otherwise everything seems correct.

Verify that $\text{Dec}\,=\,90^\circ$ is along the Z axis:

$\text{RA}\,=\,90^\circ$ is along Y.

Let's go the other way:

Two things to note: we're using the arctan2 function, which produces $\theta$ with the correct sign. We then handle negative angles specially to ensure $\text{RA} \in [0^\circ, 360^\circ]$. Check a few cases:

Convert the whole footprint:

And plot it in 3D: